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ABSTRACT 


Omega is a worldwide VLF navigation system which utilizes phase 
differencing techniques. The VLF signals of Omega are perturbed by 
phase uncertainties caused by factors such as diurnal effects which 
are quite predictable and other factors which are not so predictable. 
This dissertation models the effective phase velocity uncertainties 
which result from the less predictable phase uncertainties as jointly 
Gaussian perturbations to the phase velocities. This approach of 
modeling VLF propagation for Omega has not previously been taken. The 
more predictable factors are incorporated into determination of the 
mean values of the phase velocities. This model allows any interfre- 
quency and interstation correlation to the phase disturbances of Omega 
to be readily taken advantage of to improve the estimation performance. 
From this model, the maximum likelihood position estimation equations 
are derived for Omega with an arbitrary number of frequencies. Monte 
Carlo simulation results are reported for four-frequency Omega as SNR 
and phase velocity uncertainty are varied. These results allow quanti- 
tative error prediction for Omega utilizing different fourth frequencies 
to be made for the first time. Several fourth frequencies which can 
be easily implemented by the Omega stations, to place less stringent 
requirements on ambiguity resolution techniques, are examined. The 
requirement for resolving the inherent lane ambiguity of Omega arises 
when application of Omega to search and rescue systems such as the 
planned worldwide Global Rescue Alarm Net (GRAN) is considered.  Four- 
frequency experimental data was processed with the maximum likelihood 
estimator and the results, which demonstrate the potential of the 


estimator, are presented. 
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I. INTRODUCTION 


The Omega navigation system is being developed with the 
intent of achieving a worldwide navigational aid. ТЕ 15 
anticipated [1], that the eight station Omega system will be 
fully operational by 1976. There are, at this writing, two 
fully operational stations radiating the prescribed 10 kw 
of power at the three present Omega frequencies; 10.2 kHz, 
11 1/3 17278101 13:02КН2-г. These stations are located in 
North Dakota and Norway. The exact locations are contained 
in Table I. The Trinidad station, which will eventually be 
replaced by the station in Liberia, is operational at re- 
duced power; and the Hawaii station is expected to be fully 
operational at full power by 1975. 

The choice and number of frequencies for the Omega sys- 
tem has been predicated on the use of Omega as a navigation 
aid by ships and aircraft which would have dead reckoning 
(DR) and alternate fix capabilities. In the recent past, 
however, the United States Navy [2] and the United States 
Air Force [3], along with several other organizations, have 
been investigating the use of Omega for applications such 
as search and rescue (SAR). The requirements of a SAR Sys- 
tem, using Omega as a position location information source, 
are much more demanding than those of navigation systems. 
In general, no a priori position information is available 


at a SAR incident location. On the basis of what information 





TABLE I 


x 
OMEGA STATION LOCATIONS 


Station Latitude Longitude 
| лата, Norway 66°25'15" М 13°09'10" Е 
хх 
Trinidad 10942'06.2" N 61°38'20.3" w 


d 
Liberia (proposed) 


Haiku, Hawaii ТГ АВЕСТО EN Шар пон ри и 
La Moure, North Dakota Им парира У 98220"08,8" Ч 
Reunion ТО оз тр у ныс SL T7 ROTERE 
Trelew, Argentina 43°03'12.5" S 65°11'27.7" ы 


Australia (proposed) 
Tsushima, Japan 30 3615. 3 М 129 20749 5 НЕ 


* 
Forestport, New York 43°26'40.9" N 75 050978" W 


; Information for this table was obtained from the 
Defense Mapping Agency, Washington, D. C. 


ax | 
The Trinidad station is expected to be replaced by 


station in Liberia, 


KEK 
The Forestport station remains only as an experi- 


mental station and is not part of the Omega system. 
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can be received at the SAR incident and relayed to rescue 


forces, the SAR position must be determined. 


A. INTRODUCTORY BACKGROUND 

The first proposal to use satellite relayed Omega infor- 
mation for SAR purposes was by Samek and Pike [4]. The 
principle ТА taken from the Omega position locating experi- 
ment (OPLE) proposed by NASA's Goddard Space Flight Center. 
The OPLE was carried out in 1967[5] utilizing the ATS-1 and 3 
satellites and demonstrated that Omega navigation signals 
could be received at a remote site and retransmitted via 
satellite to a ground center which uses these signals to 
determine the location of the retransmitter. The success of 
this experiment engendered the advanced OPLE concept.  Ad- 
vanced OPLE was to use a UHF satellite link frequency (402 MHz) 
instead of the VHF link frequency (149 MHz) of OPLE, smaller 
transceivers (platforms), random access techniques for the 
satellite link, and an ambiguity resolution technique for the 
Omega information. Omega ambiguities are described later in 
this chapter. These changes to OPLE allowed for much more 
diversified anticipated usage of the concept. The random 
access technique (many users using a single satellite at 
arbitrary times) obviated the requirement to interrogate the 
platforms and thus they could be much smaller and lighter. 
This allowed for р атн application of the concept to tasks 
such as tracking of animals for zoological purposes.  Reso- 


lution capability for the ambiguity problem would mean no 


TI 





9 
S , і 


a priori position information would be required to locate 
the transceiver. The advanced OPLE concept was not carried 
out because of fund limitations and the interest of Naval 
Air Test Center (NAVAIRTESTCEN), Patuxent River, Maryland 
in using the equipment for their GRAN progran. 

1. GRAN 

In 1969, NAVAIRTESTCEN, Patuxent Ríver, Maryland, 
proposed the Global Rescue Alarm Net (GRAN) as a worldwide 
SAR system based on the advanced OPLE concept. GRAN will 
utilize hand-held transceivers called SAR communicators 
(SARCOM's) to receíve Omega signals and relay them, shifted 
in frequency and occupying reduced spectrum, via SAR satel- 
lites (SARSAT's) to ground SAR centers (SARCEN's), Figure l. 
At the SARCOM, unit identification information is added to 
the spectrum. Іп 1970, NAVAIRTESTCEN carried out experi- 
ments utilizing low transmitted powers and demonstrated the 
feasibility of low power SARCOM's for GRAN [6]. The 100 kHz 
band of frequencies from 406.0 MHz to 406.1 MHz was assigned 
for use in a SAR satellite link such as GRAN at the World 
Administrative Radio Conference for Space Telecommunications 
(WARC-ST) in Geneva in 1971 [2]. 

It is ДЕК с уаз” by the originators that GRAN will 
serve all Soe GRAN will be available to any military 
and civilian users who purchase a 5АКСОМ, In keeping with 
this concept, the cost of the SARCOM units must be kept to 
the minimum possible level. This constraint prohibits pre- 


processing of the received Omega information in the SARCOM, 
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| 





'uo3sÁS куно PUL “Т э10814 


NOIIVIS VIHHO 4 
Та Сан LINN ROOWVS 


NHOWVS 


N 
хаттап N 
NOIIVIS VIMO | A 
xur[uaoq N 


ad LVSUVS 


ини! 


13 





In order to meet the portability requirements necessary so 
pleasure craft and perhaps even hikers can carry the units, 
weight and size of the SARCOM units must also be kept to a 
minimum. Studies by NASA [7] have demonstrated the feasi- 


bility of the following characteristics for the SARCOM's: 


Volume 1000 i 

Weight l kg 

Output Power 5 watts maximum 

Frequency 406.0 - 406.1 MHz 

Bandwidth 249 kHZ 

состоя Period 3 min. broadcast, repeated as 


battery power allows 

Identification ¡Data 36 (bits) - sufficient for 
nine digit number such as 
a social security number 


The SARCOM units will operate sequentially in three 


modes when activated. The modes are: 


Mode 1. Acquisition 
Mode 2, Identification Data 
Mode 3. Omega Data 


In mode one, the SARCOM will transmit only an internally 
generated Acquisition/Recognition (A/R) tone in the allotted 
band of frequencies for a period of about 12 seconds. The 
maximum power capability of the unit will be concentrated 

in this tone to allow the SARCEN phase-locked-loop receivers 
to acquire the retransmitted tone and switch to a narrow 


tracking bandwidth for signal-to-noise ratio improvement. 
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The second mode is then commenced automatically and SARCOM 
unit identification information is transmitted digitally by 
phase-shift-keying the A/R tone with a code preassigned to 
the SARCOM unit. The second mode will consist of about six 
seconds of epp ES In the final mode, the received 
Omega information is retransmitted along with the A/R tone. 
The total transmitted power of the SARCOM will be divided 
between the A/R tone (approximately .5 w) and the retrans- 
mitted Omega signals (remaining power equally divided). 

This third mode is expected to last for approximately 180 
seconds, which would allow for 18 of the ten second Omega 
format periods (see the following section) to be relayed. 

A possible spectrum of signals to be retransmitted if four- 
frequency Omega is adopted and the fourth signal is at 10.88 
kHz is shown in Figure 2. The signals in the spectrum cor- 
respond to frequency shifted versions of the Omega signals; 
signal number one corresponds to the 10.2 kHz signal, signal 
Ето to the 10.88 kHz signal, signal three to the 11 1/3 kHz 
signal and signal four to the 13.6 kHz signal, The actual 
number of different frequency Omega signals and the options 
available in choosing these frequencies is discussed in de- 
tail later in this thesis. The center frequency of the 
spectrum in Figure 2 is referred to only as fo. This is done 
since to utilize the entire 100 kHz bandwidth of the SARSAT 
as well as to avoid concurrent identical spectra to be trans- 
mitted from several SARCOM's, the 2.5 kHz SARCOM spectra 


will have to be spread throughout the 406.0 MHz - 406.1 MHz 
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band of the SARSAT. This may be accomplished by assigning 
different frequency A/R tones to different users. Reference 
[8] contains additional discussion in this area. The spac- 
ing of the retransmitted Omega signals relative to the А/В 
tone will be accurately controlled by the SARCOM. Thus, 
when the phase-locked-loop receiver in the SARCEN locks on 
and tracks the A/R tone, the Omega signals may be recovered 
with essentially all of the frequency drift and doppler 
shifts which affect the spectrum, after its creation in the 
SARCOM, eliminated. 

The SARSAT system to be utilized for GRAN will con- 
sist of a system of three to four geosynchronous satellites 
and additional polar orbiting satellites. The exact number 
and types of orbits of the satellites is still being con- 
sidered by NAVAIRTESTCEN. A discussion of the requirements 
for the satellite system as well as link calculations is 
contained in [8]. 

The SARCEN system will consist of several SARCEN's 
in order to allow some SARCEN to continuously be in the 
coverage of 2001 SARSAT. Each SARCEN will search the pre- 
assigned band of frequencies with a set of phase-locked-loop 
(PLL) receivers. The bandwidth of the PLL in the search 
mode will be larger than when tracking a SARCOM relayed 
message in order for the acquisition time to be small. 

When an A/R tone has been acquired, the bandwidth of the PLL 


will be decreased to improve the SNR from the PLL and allow 
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lock to be maintained after the power of the A/R tone is 
decreased for the third mode of SARCOM transmission. 

After the SARCEN has acquired a retransmitted SAR 
message, identified the sending SARCOM from the identifica- 
tion data, and recorded the retransmitted Omega signals, 
the position of the SAR incident must be determined. The 
optimum method for locating the retransmitting SARCOM posi- 
tion is determined by this thesis. 

2. Omega 
а. History of Omega 

Following World War II, the MIT Radiation Labora- 
tory which had been actively developing military navigation 
aids, was closed [9]. Loran A and Loran C, as presently 
used, were results of their work. A small group from Divi- 
sion 11 of the Laboratory moved to Cruft Laboratory, Harvard 
University. In 1947, Professor J. A. Pierce of this group 
proposed Radux [4] as a navigation system. The technique 
of Radux was very similar to Loran except that the trans- 
mitted signals were phase modulated with a 200 Hz signal. 
Because of SEC 72201, stable phase characteristics of 
signals in the proposed frequency range of 40 kHz to 50 kHz, 
it was expected that the modulation could be used in phase 
comparison with quite accurate results. 

The results of experiments with Radux by the 
Navy Electronics Laboratory indicated that the accuracy of 
Radux could resolve ambiguities which would result from a 


10 kHz signal. This, along with the fact that VLF 
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frequencies were relatively unusable for communications pur- 
poses, resulted in a composite system Radux-Omega which 
radiated the 40 kHz signals of Radux interrupted by 10 kHz 
bursts of a coherent signal. The 10 kHz signal yielded 
nominally eight nmi (nautical miles) ambiguities which could 
be resolved by the 200 Hz modulation of the 40 kHz signal. 

It was soon found that this system was not as 
good as could be obtained. The 10 kHz signal had a much 
longer usable range than the 40 kHz signal and could thus 
have much longer baselines (the great circle connecting two 
Omega stations) thus improving the accuracy of the system 
markedly. 

b. Present Omega 

Several different frequencies in the 10-14 kHz 
range were used in experiments with Omega. The frequencies 
decided upon for final implementation were 10.2 kHz, 11 1/3 
kHz, and 13.6 kHz common to all Omega stations and other 
frequencies unique to the stations to be transmitted for 
other than phase difference navigation. A total of eight 
stations TE decided upon to achieve worldvide coverage 
with a station redundancy sufficient to assure an adequate 
fix anywhere in the world. The three common frequencies are 
radiated as on-off-keyed, time-division-multiplexed signals 
in a ten second format as shown in Figure 3. The nominally 
one-second pulses are all separated by 0.2 seconds. Each 
station transmits the common frequencies in the same sequence, 


but starting at a different time. 
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Omega, being a phase differencing navigation sys- 
tem is necessarily a hyperbolic system. That is, the isophase 
contours of two stations (lines on the earth for which the 
phase difference of two common frequency signals from two 
stations is constant) form a set of hyperbolas shown in 
Figure 4, which has as foci the two stations radiating the 
signals. As shown in Figure 4, two adjacent hyperbolic lines, 
henceforth called "lines of position" or LOP, separated by a 
distance d, on the baseline can be shown to be spaced by 


1 
d..= d, csc (0/2) off the baseline [10]. The angle O is that 


27 1 
angle formed by the two great circles passing through the 
stations and a point midway between the adjacent LOP off the 
baseline. The spacing of two LOP off the baseline of two 
stations then will be equal to the spacing of these LOP on 
the baseline multiplied by csc (8/2). 

Because it is a phase difference system, Omega 
has naturally occurring ambiguities. Оп the baseline of 
two stations, the sane phase difference of a common frequency 
f radiated from the stations will be seen every 23 where À 
is the wavelength of the frequency f. This results in lanes 
between unambiguous LOP of width Е: A on the baseline. Off 


2 
the baseline this spacing will be + A ese (0/2). 
This concept is easily extended to multiple fre- 
quency transmissions by considering a signal of frequency - 
which is a common clock frequency to all Omega stations and 


is utilized to create the K transmitted signals - (17:27:01 17:24 


e è e К). In the present Omega system K is three. In order 
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for phase comparison of the common frequencies from the sta- 
mons to Бе useful in locating a remote receiver site, the 
frequency f. at all stations are maintained very nearly 
phase coherent. The ith Omega signal to be transmitted is 
created from division of Е. by an integer ni. The baseline 
lane width associated with frequency г. will be = vs/ 2f, 
where Е 8 is the largest common multiple frequency of the 
transmitted frequencies. This is the smallest frequency 
which can be obtained by frequency differencing techniques 
from the K transmitted frequency. The value of Mes the 
effective phase velocity of Е o? is determined by the phase 
velocities of the K signals. However, Ns is nominally equal 
to c, the velocity of light in vacuum, so we can get a nomi- 
nal value for Sh. by assuming MS Š c. The term "unambiguous 


lane will be used in thils® thesis to refer to 3) The 


ар 
nominal baseline size of the unambiguous lane may be deter- 
mined by the folloving procedure. The integer divisors, 

na (i = l1, 2, sse, K) are first separated into the, products 
of powers of their prime number components (see example which 
follows). The product of the largest powers of all prime 
numbers which occur, then gíves the divisor п. which yields 
f from a f e/n.. Baseline unambiguous lane width is then 
~ c/2f. = с • n /2f . The present value of f (816 kHz) used 


by the Omega stations allows a simple rule of thumb to find 


a nominal baseline unambiguous lane vidth easily from nc. 


161,875.2(nmi/sec)* M 





l Е ~ Р ~ 
5^, (nni) =,0992 D. по /10 (1) 


2 + 816,000 Hz 
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The present three frequencies transmitted by the Omega sys- 
Bew; 10.2 kHz, 11 1/3 ЕН2, and 13.6 kHz yield n,- iS 
44720 в1псе п 5 287) 5 180, по = 23.32%4 72, апа DNE 
22. 3 «5 = 60, Thus from (1) the baseline width of the 
unambiguous lanes is nominally 72 nmi. 
3. VLF Propagation 

Since the Omega system depends upon phase comparison, 
unexpected phase shifts in the propagation path cause posi- 
tion estimation errors. The causes of the anomalous phase 
shifts must therefore be considered. For the Omega systen, 
frequencies in the 10-14 kHz region have been chosen because 
of the extremely long ranges which VLF signals propagate and 
because of the remarkable phase stability of these signals. 
The, propagation problem at VLF has been studied both experi- 
mentally and theoretically by many researchers for many years. 
References [11, 12, 13, and 14] contain very detailed and 
thorough discussions of VLF propagation. 

In general, for distances less than several Megameters 
(Mm) from the transmitter, the easiest theoretical analysis 
of VLF propagation can be obtained by a combination of ground 
and sky wave contributions. The ground wave provides the 
field which would be seen in the absence of the ionosphere 
and its attenuation with distance from the transmitter is 
generally attributable Во En curvature [11]. The sky wave 
is the transmitted wave which is reflected in the earth iono- 
sphere cavity. The contribution to the total field from 


skywaves can be determined by summing all of the skywaves 





present at any distance from the transmitter using ray 
theory [11]. At distances beyond a few Mm, the number of 
skywave contributions becomes unwieldy and waveguide theory 
of propagation becomes much easier to use. 

Waveguide theory for VLF is generally used for mak- 
ing theoretical predictions for Omega. Propagation predic- 
tions (to be discussed later) do attempt to predict phases 
of signals at ranges where ground wave effects are signifi- 
cant; however, Omega is considered most predictable and 
thus most reliable at distances in excess of 2-3 Mm from 
the transmitters. At these distances, single mode waveguide 
theory seems to describe the propagation quite well. 

Waveguide theory of propagation describes the earth 
and ionosphere as forming a spherical waveguide. The height 
of the ionosphere above the earth, nominally 90 km at night 
and 70 km during the day, is on the order of several wave- 
lengths for the 10-14 kHz signals of Omega and thus wave- 
guide theory would seem to apply if the boundaries of the 
waveguide could be adequately described. Тї 15 generally 
accepted [15], that the propagation effects causing phase 
variation, which is what is of concern for a phase comparison 
System such as Omega, can be divided into two categories, 
temporal and spatial. Among the temporal variations, the 
diurnal variations have the most dramatic affect and are 
probably one of the most predictable factors. This variation 
is in general accounted for by using different effective 


phase velocities for each frequency for dark and for lighted 
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transmission paths. Professor J. A. Pierce of Harvard Uni- 
versity has conjectured that adequate position estimation 
using Omega may be attainable using the single waveguide 
node phase vetoeities given in Pable IT. These values were 
obtained from [16] and from personal communication with 
Professor Pierce. Other temporal variations are long term 
dependence on solar cycle, lunar, and annual variations. 

The spatial variations of importance appear to vary 
with propagation bearing, earth conductivity over the trans- 
mission path, geomagnetic latitude, and the ionosphere con- 
ditions. The magnetic field has an apparent marked effect 
on the ionosphere and thus creates phase variations depen- 
dent upon the transmission path. These spatial variation 
causes are not independent quantities and thus must be 
treated together [15]. 

The Defense Mapping Agency attempts to account for 
all of the above phase variation sources as well as other 
lesser known sources in the Omega Propagation Correction 
Tables which are published regularly for the 10.2 kHz and 
13.6 kHz Omega signals. The computer programs which accom- 
plish this are described in [15]. Additionally, force 
fitting of the corrections to actual data is continually 
carried out. The use of these tables is simple; knowing the 
date and the approximate receiver location, phase correction 
factors to be added to the observed phase of a known fre- 
quency signal from a known station are looked up in pre- 


published tables. After these corrections are made, the 
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f (kHz) 


10.2 
11 1/3 
13.6 


10.88 


TABLE II 


ВЕРОМ РЕНЕ ОГ MODE PHASE VELOCITY RATIOS 


vic Day 


1.00270 
1.00139 
9 906 5 


12.009187 
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v/c Night 


„99960 
„99869 
99731 


„99902 





phas difference for a particular station pair and frequency 
choice is computed and the corresponding LOP is determined 
from specially prepared Omega charts. Figure 5 shows a 
sample section of an Omega chart. These charts are prepared 
based upon a phase velocity 1.00261 times the speed of light 
Mmveacuum and thus the phase corrections attempt to correct 
the phase to this effective phase velocity in determining 


the receiver location. 


B. PROBLEM STATEMENT AND APPROACH 

This dissertation determines the optimum procedure to 
establish an Oméga receiver's location. In order to accom- 
plish this task, a model for the Omega problem was first 
formulated in which the effective phase velocities of the 
transmitted Omega signals were assumed known. Since true 
phase velocities of the signals are continuously varying 
along the transmission path as ionospheric densities and 
many other factors change, effective phase velocity is used. 


This velocity is that which satisfies the equation 


у = —— (2) 


where d is the distance measured on the globe which the 
signal travels, w is the angular frequency of the signal, 
and $ is the observed phase in radians of the signal at a 
distance d trom the transmitter. 

Utilizing this assumption and the assumptions about the 
noise which will be given in Chapter II, the maximum likeli- 


hood LOP estimation equation for two-station, multiple 
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frequency Omega was derived. This derivation is presented 
mU ppandix C. Using this result and the effective phase 
velocities for single mode propagation given in Table II, 
some four-frequency experimental Omega data was processed. 
The data was collected METRE to the four frequency 
Omega experiment described in detail in Chapter VI. The re- 
sults of this effort were no lane errors (defined in Chapter 
I for transmission paths in daylight, but a significant 
number of lane errors for dark or partially dark transmis- 
sion paths. The conclusion which was reached then was that 
the Table II effective phase velocities were not sufficient 
predictions of phase velocities actually being experienced. 
Evounting for diurnal effects only, in predicting effective 
Ese velocities, and assuming single mode propagation was 
Doc apparently justified. | 

A further search into the literature revealed that at 
reasonably long ranges from VLF transmitters, predicted 
phase velocities using propagation prediction methods dis- 
agreed with experimental results by a few parts in 105171 
This indicated that the best available methods of predicting 
the disturbances to the phase of VLF signals for typical 
distances could yield errors in phase differences of signals 
. from two stations on the order of one radian or more. This 
Beale zation Med то the model which was used in this’dissier- 
bation. 

Since phase is related to effective phase velocity by 


(2), phase disturbances can equivalently be considered as 
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due to effective phase welocity uncertainties. The many 
effects which cause phase uncertainties and thus effective 
phase velocity uncertainties imply, as a result of the cen- 
tral limit theorem, a model of the phase velocities of 
jointly Gaussian random variables. This new approach not 
previously considered in the literature allows any correla- 
tion between the phase disturbances of signals from the 
same station or different stations to be taken into account 
in the estimation process. The effective phase velocity 
can be predicted as accurately as necessary or possible for 
the different frequency signals and these values used for 
the mean values of the jointly Gaussian phase velocities. 
The entries in the covariance matrix of the phase velocities 
reflect the uncertainties of these predictions and the 
correlation of the uncertainties. This is the first model 
for position estimation which allows advantage to be taken 
of the correlation between different frequency signals from 
an Omega station. 

References [18, 19, and 20] reveal a significant correla- 
tion between phase perturbations of Omega signals transmitted 
from a single station. The interstation correlation is in 
general considered to be negligible [19]. If some appro- 
priate standard deviations of the phase velocities on the 
order of a few parts in ТО ате assumed and the correlation 
between all signals are assumed known, the probability den- 
sity functions of the assumed jointly Gaussian phase velo- 


cities are completely known. These parameters must be 
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ши mentally determined. Professor J. A. Pierce of Harvard 
University is presently investigating the correlation of 
phase disturbances of Omega signals. 

Utilizing this new approach and the assumption of addi- 
те vhite noise, which is justified in Chapter II, the maxi- 
mum likelihood estimation equations for Omega were derived 
(Chapter IV) for high SNR. The performance of the resultant 
estimator for various fourth frequencies which could be 
easily added to the present three-frequency Omega format was 
investigated by Monte Carlo simulation. The results of these 
Sumi lations (seel Chapter V) allow for the first time, a 
quantitative measure, in terms of expected errors, of the 
performance of Omega with different fourth frequencies which 
could be added to Omega to yield larger unambiguous lanes. 
The maximum likelihood estimator was used to process experi- 
mental four-frequency data (Chapter VI) and the performance 
was found to conform in general to the simulation results, 
Preliminary comparison of the estimator performance indicates 
it is superior in terms of lane errors (defined in Chapter 
IV), the most important consideration for SAR systems parti- 
cularly, to any other known estimation procedure. 

The equation for the maximum likelihood estimator re- 
quires prior estimates of the doppler shifts, signal ampli- 
tudes and noise spectral density heights for the Omega 
signals received by the SARCOM. The maximum likelihood esti- 
mation equations for these quantities are contained in 


Appendices A and B. The results of application of the signal 
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amplitude and noise density estimation equations to experi- 
mental data are presented in Appendix B also. These equa- 
Brons allow optimum estimates of signal-to-noise ratios to 
Ba found for Omega signals. The optimum method for determin- 
ing SNR of Omega signals di not been presented before. 

A four frequency experiment utilizing a satellite link 
is planned by Naval Air Test Center, Patuxent River, Maryland 
for late 1974. A prototype, portable transceiver is being 
constructed and will be transported to various sites in 
Europe and the eastern United States for relay of four fre- 
quency Omega signals. The LES-6 satellite will be used. 
The processor for the ground station which has been selected 
for use is the maximum likelihood estimation program pre- 


sented in Appendix E of this dissertation. 


CE ORGANIZATION OF THESIS 

Chapter II of this thesis presents the mathematical 
model of the estimation problem. The assumptions necessary 
to describe the problem so that a solution can be obtained 
and their justification are presented. 

Chapter III addresses the resolution of the ambiguity 
problem which arises in a differencing system. The size of 
the unambiguous lanes, determined by the number of signals 
and the choice of their frequencies, determines the require- 
ments from an ambiguity resolution technique. Several 
potential techniques and their anticipated capabilities, 
EUch wibstnus determine frequency choice for a SAR system, 


are presented. 
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Chapter IV contains the derivation of the maximum likeli- 
hood estímate of the LOP resulting from two stations each 
transmitting an arbitrary number, K, of common frequencies. 
Since there are only eight locations in the Omega format for 
the nominally one second pulses of different frequencies, K 
will necessarily be less than or equal to eight unless more 
than one frequency is allowed to occupy a pulse by time shar- 
ing. The resultant estimation equation when K is set equal 
to four, the number of frequencies presently e considered 
by. the managers of GRAN, is presented. The estimation equa- 
tion which resulted when deterministic, known phase velocities 
were assumed is addressed in Appendix C. 

Chapter V presents the results of Monte Carlo simulation 
which was carried out for several easily implemented fourth 
frequency choices, an expected range of variances of the 
phase velocities, and a range of SNR. Тһе simulations were 
carried out to aid in selection of a fourth frequency for 
Omega which will allow an ambiguity resolution technique to 
resolve the unambiguous lane problem. The trade off in per- 
formance in terms of estimation error which is experienced 
as unambiguous lane size is increased as a function of the 
necessary parameters is presented. Although the addition of 
one frequency to the Omega format is expected to be accept- 
able to the program managers as well as the Omega station 
host countries, additional frequencies and/or modulation of 
a frequency has been stated to be unacceptable by the Omega 


program office and is, therefore, not considered. 
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Chapter VI presents the results of processing four-fre- 
quency Omega data from two stations with the maximum likeli- 
hood estimator. The data for the experiment was collected 
by Texas Instruments Corporation of Dallas, Texas for 
NAVAIRTESTCEN in the fall of 1973. The fourth frequency 
which was transmitted by the two experimental Omega stations 
Sirinidad and Forestport, New York) in addition to the 
standard three Omega signal frequencies was 10.88 kHz. A 
discussion of the results and an explanation of the errors 
which were experienced is presented. The estimation was 
done first by assuming the effective phase velocities given 
in Table II for single mode propagation and second by using 
propagation correction factors provided by the Defense Map- 
ping Agency. A comparison of the results is presented. 

Chapter VII presents the derivation of the multiple sta- 
tion estimation equation and the prospects of being able to 
implement its computer solution. An alternative approach 
to the multiple station problem which requires less computa- 
tion time but is nonoptimum is presented. 

Chapter VIII presents the conclusions drawn from the pre- 
vious chapters. Several appendices contain important methods 
and results which are not required in the general flow of 


the chapters. 
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Шинэ STHEMATTEAT MODEL 


Consider two Omega stations each transmitting K frequen- 
cies occupying K of T nominally one second positions 
in the predetermined format shown in Figure 3. Let Y(t) be 
a 2K dimensional column vector, where the first K components 
represent the signals received from station one and the last 
K components represent the signals received from station two. 


We shall represent Y(t) as 
мге 50721150: (3) 


where S(t) represents the signal portion of the received vec- 
tor and N(t) represents an additive noise process. We can 


express the ith component of S(t) in (3) as 
si(t) = A, cos ACA > А ЈЕ - 9, + 9.1 (4) 


where А. isi the апр] иде о the ith signal, wj is the angu- 
Шат frequency of the itti signal, Az is the doppler frequency 
shift of the ith signal, 9, is the propagation phase delay 

ЕЙ Ере ith signal and 9. is the contribution of other than 
propagation delay to the phase of the ith signal. The com- 
ponents of the 2K dimensional signal vector will be ordered 
so that the angular frequency wj equals the angular frequency 
Since quite accurate common frequency signals from the 


ICE" 


two Omega stations are necessary for navigation, кі = W; 


will be assumed. 
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For GRAN, an acquisition and reference signal (A/R tone) 
will be generated in the SARCOM and transmitted with the 
Omega data. The transmitted information will occupy 2.5 kHz, 
Figure 2, a smaller bandwidth than the 3.4 kHz received 
signal bandwidth, for spectrum conservation and will be 
centered at a frequency in the allocated band of frequencies 
from 406.0 - 406.1 MHz. The A/R tone will serve two func- 
tions. For a short period at the beginning of the transmis- 
sion period it will be radiated alone by the SARCOM so the 
SARCEN phase-locked-loop receiver can achieve lock and track 
the subsequent transmission. The other function of the A/R 
pone” wus) to allow phase shifts and frequency drifts in the 
transmission channel through the SARSAT to the SARCEN to 
pe eliminated. The particular contributors to these effects 
will be satellite link delays and frequency drift in the 
satellite frequency translator. Proper processing at the 
SARCEN will then result in A. values in (4) being caused 
only by relative motion of the SARCOM with respect to the 
Omega stations. The frequency drifts in the SARCEN or 
Omega stations during the transmission period will be 
negligible because of the accurate frequency standards which 
will be used. 

The same transmission channel will be used to receive, 
relay, and process the (Kti)th frequency as the ith since 
they must be at the same frequency. Therefore, with the 


exception of timing errors at the stations, the phase shift 





Ода = 0 1» илэг The phase shifts in the К 
different channels will result from different phase shifts 
at different frequencies in the transmission channel from 
the SARCOM to the SARCEN and from the frequency division in 
the SARCOM. The compacting of the spectrum in the SARCOM 
must be accomplished by multiplication in each receiving 
channel with an appropriate frequency signal created by fre- 
Buency divisfon of a local oscillator signal. The frequency 
division will be done by digital counting circuits whose 
initial counting level will be a random quantity. The phase 
shifts in the K channels which result will thus be mutually 
independent random variables, uniformly distributed on the 
interval (-T, T]. It can be shown that this assumption is 
morsaffeeted by the addition of constant phase shifts to the 
signals due to Omega station timing errors since a uniformly 
distributed random variable plus a constant remains uniformly 
distributed. 

The propagation phase shifts of the signals from station 


one and station two will be written as 


d 
$. = m, vy cm, en) 
O 
and 
$ „|| an 27m ШШ cx CSD 
к-і 1 Vary 002 К+1 = = : 
O 


an (5), У1( ) is the effective phase velocity of the ith 


signal from the Omega station to the SARCOM receiver, These 
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phase velocities are written as a function of the distance 
from the transmitter to allow for incorporation of changes 
in the predicted value of the phase velocities, e.g., from 
daylight path to dark path conditions. The integers m, are 
the number of whole ОО of phase delay the ith signal 
has experienced from the Omega transmitters to the SARCOM 
receiver. The values x] and x, are the great circle dis- 
tances from the respective Omega stations to the SARCOM. 

It is convenient since only an LOP estimate can result 
from two stations, to estimate only one variable rather than 
two (хү and хо) in (5). Thefgraphical description of this 
method is contained in Chapter IV. It consists of defining 


вое 5 чм хх + ох, = І and replacing? x. by x. Then, by holding 


T 2 


the value of L constant the only variable to estimate is x. 


l 


Substitutions) of these definitions into (5) yield 


x 
ж So АИ s : 
$i m | У (ү) 21m, 1 < :1< Кк (а) 
о 
and 
JE 
d 
фул. = б, ? 一 一 二 一 一 -2Tm XL JI 07271620. 
K+i i Урд Ч ү) K+i 
x 


Equations (5) and (6) are written as two equations in order 
to emphasize that We = Opis and to allow the limits of the 


integrals to be written. The value of x is the great circle 


distance from the first station to the receiver at the 
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beginning of the data collection period (t-0) and L is the 
Sumeof thee great circle distances to the two stations for 

t = 0. Any SARCOM motion must be accounted for by using the 
doppler estimation information from Appendix A and the posi- 
tion estimate. | 

The large range over which VLF signal strengths can be 
expected to vary and the impulsive nature of the atmospheric 
noise at VLF frequencies [11] dictate the use of a limiter 
in the receiver. This thesis will assume that the receiver 
will use a hard limiter which will yield the necessary 
dynamic range as well as minimize the effect of the impulsive 
noise, 

The components of the noise column vector N(t) in (3) 
will be assumed to be sample functions from white Gaussian 
random processes, the ith component of N(t) having two-sided 
noise spectral density height of N,/2. The ith and (K-ti)th 
noise components will be received in the same channel of the 
receiver Since they are associated with signals of equal 
frequency. But,as shown by Figure 3,they are received at 
different times and will, therefore be uncorrelated. The 
polse willfactually be thle sum of thermal and atmospheric 
noise and the rationale for the Gaussian assumption is as 
follows. The atmospheric noise, expected to dominate at the 
VLF frequencies, will be the sum of many atmospherics from 
around the globe. The distant noise impulses will have in- 


sufficient amplitude to dominate the hard limiter in the 


receiver, but nearby impulses will dominate the hard limiter. 
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Since VLF signals propagate very well over long distances, 
the small amplitude impulses of noise will occur quite fre- 
quently and by the central limit theorem their sum will 
approach a Gaussian random process. The strong impulses 
which dominate the hard limiter, yielding approximately zero 
SNR to the receiver during their existence, will occur much 
less frequently and with a wise choice of bandwidth for the 
K channels of the receiver, each impulse will capture the 
hard limiter for only a small time. It can be shown that 
this time will be on the order of the inverse of the band- 
msedth Of the filter, so that a 100 Hz filter will allow the 
hard limiter to be captured for a time of about 0.01 sec. 
by an impulse. 

The 2K dimensional phase velocity vector V will be used 
to represent the effective phase velocities which are 
modeled as jointly Gaussian random variables. The covariance 
matrix of the components of V will be written as K Since 
inadequate evidence is presently available to justify differ- 
ent uncertainties for the effective phase velocities when 
modeled in this fashion, it is convenient for the simulations 


reported on in Chapter V to assume 


р 0 
Т 2 
ККЕ о 2107) 
0 р 


Ша (7), 9 * -is the scalar variance, assumed equal for all 


signals. The symbol p represents the K-by-K correlation 





matrix for the signals from a single station. The symbol O 
represents a K-by-K zero matrix, indicating interstation 
correlation is negligible. The rationale and justification 


for these assumptions is contained in Chapter I. 
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DPI AMBIGUITY RESOLUTION 


Phase differencing systems such as Omega necessarily are 
ambiguous. The size of E. unambiguous lanes (the lanes 
formed by indistinguishable LOP from a station pair) is deter- 
mined by the choice of frequencies. A discussion of how the 
frequencies determine the lane size is contained in Chapter Г. 

The Radux navigation system 19) which immediately pre- 
ceded the Omega system had a 200 Hz modulation of its carrier 
frequency, 40 kHz. This allowed for nominally 400 nmi un- 
M ous laná which were considered to be adequate for navi- 
gation. It was the short baselines, due to the short 
propagation ranges, which caused Radux to be abandoned and 
not its ambiguities. The short baselines caused the diver- 
penccror the ГОР off the baseline to occur very rapidiy as 
distances from the baseline increased, thus reducing the 
accuracy of the system, When Omega was in the design stage, 
dead reckoning (DR) systems and the reliability of receivers 
were improving markedly. This influenced the final design 
Ша Опера to concentrate) on one ог two frequencies which 
would be used for navigation. Many Omega navigators presently 
use only the 10.2 kHz Omega signal. Also, the Defense Mapping 
Agency only publishes propagation prediction correction 
factors for the 10.2 kHz and 13.6 kHz Omega signals. As 
shown in Chapter I, the three navigation signals presently 
radiated by the Omega stations yield nominally 72 nmi un- 


ambiguous lanes at the baseline. 
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This has been adequate for all systems using Omega to 
date. The SAR applications such as GRAN, which prompted 
this study, however, require an ambiguity resolution tech- 
вое in order to be successful. Since the retransmitter 
at the SAR incident will in general have no prior knowledge 
of position which can be relayed to the SAR forces, the 
assumption must be made that complete position estimation 
must be accomplished from the retransmitted Omega data. 
There are presently five methods which are seriously being 
investigated as potential ways to solve the ambiguity 
problem. They are: 1) addition of frequencies to the Omega 
format, 2) pulse time-of-arrival measurement, 3) ratio of 
signal-to-noise ratios, i.e., signal-to-signal ratios, 4) 
multiple LOP examination, and 5) modulation of one of the 
Omega signals with a large period tone. The fifth technique 
has been stated by the Omega Project Office in Naval Elec- 
tronic Systems Command as being unacceptable to the host 
countries for the Omega stations and will not be considered 
here in detail. The United States Air Force, however, is 
funding an experiment which Cincinnati Electronics Corpora- 
tion and NELC of San Diego are to carry out by 1975 which 
will modulate one of the signals from the experimental 
Emera station in Forestpoit. New York with a 226 2/3 Hz tone, 
Since the Omega Project Office has stated that a final imple- 
mentation using modulation on an Omega signal is not feasible, 


the experiment appears to be academic. 
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As discussed in Chapter I, the choice of frequencies 
determines the size of the unambiguous lanes. It is reason- 
ШЕГЕРЕ теп to expect thdadtwdifferent or additional frequency 
signals for Omega could yield larger unambiguous lanes. It 
is very likely, however, that the three present frequencies 
must remain and that, for political reasons, at most only 
one more frequency common to all stations can be added to 
the Omega format. This statement is made because of the 
many existing receivers designed for the three existing 
Omega frequencies and the great demand for the remaining 
Omega format slots for other applications. Chapter V of 
this eser а. is dedicated to the study of the perform- 
ance of four-frequency Omega for several easily implemented 
fourth frequencies. One additional frequency is not ex- 
pected to be able to completely solve the ambiguity problem; 
however, it can place much less demand on other ambiguity 
resolution techniques by expanding the unambiguous lane 
size. 

Reference [21] discusses the use of the relative times- 
of-arrival of the Omega pulses from station pairs at a re- 
ceiver to determine a lane of LOP which contains the receiver. 
One of the conclusions of [21] is that the existing three- 
frequency Omega, 72 nmi, baseline ambiguity cannot be re- 
solved in high noise areas with .95 confidence with less than 
on the order of one hour of data. The pulse time-of-arrival 
technique is, however, expected to be able to resolve nominal 


360 nmi baseline lanes with three minutes of Omega data. 
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These Tanes result from four-frequency Omega utilízing 10.88 
TA addition to the existing Omega tones. This resolution 
is expected even for worst case noise, 

Reterence |22] reports on tan experiment which examined 
the potential Omega lane size which could be resolved using 
the relative powers of the signals received from a station 
pair. The OPLE ground station in Dallas, Texas was utilized 
in measuring satellite relayed Omega data retransmitted by a 
bread-boarded SARCOM from various locations on the globe. 

It was assumed that the noise averaged to a constant over 
the essentially concurrent measurement interval for signals 
from two stations, The ground station was used to estimate 
the SNR from the two stations and with the above assumption 
eoul the noise, the ratio of the SNR's from the stations 
for a given frequency is a signal power ratio. The conclu- 
sion of [22] was that signal-to-signal power ratios E 
feasible for use in ambiguity resolution, but the limited 
experiment performed did not allow for a conclusion about 
the lane size which could be resolved. The SNR for the ex- 
periment was measured at the ground station by a non-optimal 
technique and it is possible that the maximum likelihood 
estimate of SNR derived in Appendix B of this dissertation 
could be used to advantage if this technique is applied in a 
final system. 1 

Reference [23] reports on simulations which were carried 
but Lor three - frequency multiple station Omega utilizing 


multiple LOP intersections to resolve ambiguity. А11 the 





polygons formed by LOP from stations treated pairwise are 
examined for an area of radius 250 nmi on the globe and the 
one which maximizes a likelihood function is selected as 
being the most likely position estimate. The likelihood 
function is formed by assuming additive Gaussian phase un- 
certainties affecting the phase differences which correspond 
to the LOP estimates of the receiver location when the Omega 
stations are treated pairwise. The conclusion drawn from 
the simulations is that the multiple LOP method of ambiguity 
resolution does have promise in terms of resolving the Omega 
ambiguities if multiple position estimates along with their 
associated relative probabilities are taken as a solution 
set. 

Additional information concerning the simulations reported 
Миро ОЭ is contained m [245, 25]. Lt can be fottnd fron 
these references that the standard deviations used for the 
additive Gaussian phase uncertainties for the simulations 
was 0.003291 rad for three stations and 0.010408 rad for the 
remainder, This corresponds to LOP estimate standard devia- 
tions on the order of 0.0047 to 0.0142. These values are 
said to be derived from a 10 db-Hz SNR expected from three 
stations and O db-Hz SNR from the remaining stations. Те 
was found by experiment and is discussed in a later chapter 
of this thesis that standard deviations of the LOP estimate 
on the order of 0.1 to 0.2 baseline nmi are experienced with 


daylight paths and SNR in excess of 10 db-Hz for a station 





puero using the optimum estimator. The equation 


ET = c (21/360) (8) 


allows these standard deviations to be converted to standard 
deviation of phase ЕЕЕ а 2т radians of phase dif- 
ference corresponds to one unambiguous lane on the baseline 
(nominally 360 nmi for the four frequencies which were used 
in the experiment). Therefore, standard deviations of phase 
difference uncertainties which were experienced in the experi- 
ment were, for best conditions, on the order of 0.002 to 
0.004 radians. When unambiguous lane size is increased by 
the addition of a fourth frequency as was the case for the 
experiment, the same number of polygons as were considered 
in [23] occur in a larger area on the globe. If this area 
is expanded sufficiently both by addition of frequencies and 
consideration of more polygons, the multiple LOP method of 
ambiguity resolution may be a viable technique for SAR sys- 
tens. 

The pulse time-of-arrival technique and the multiple LOP 
technique ner: to have the most promise as ambiguity reso- 
lution techniques which can be implemented. The multiple 
LOP is an alternative, although nonoptimal, technique for 
treating the multiple station problem to be discussed in 
Chapter VII. The conservation of computer time which it 
. allows would be a reason for adopting it as a method over 
the optimum estimator. An estimation procedure using the 


multiple LOP technique would be to obtain the two-station 





LOP estimates optimally and determine the set of most likely 
position estimates from a multiple LOP approach as described 
in [23]. Use of pulse time-of-arrival information with this 
approach will be very helpful in reducing the area on the 


globe over which LOP intersections must be examined. 





ТУ. TWO-STATION LOP ESTIMATION 


А.  TWO-STATION ESTIMATION FOR K FREQUENCIES 

In order to establish an estimation problem which will 
yield an unambiguous estimate of LOP for two stations and 
K frequencies, we assume that an ambiguity resolution tech- 
nique has resolved the lane ambiguity. Thus, the estimation 
is performed within one unambiguous lane. Within the pre- 
supposed unambiguous lane we can only estimate a LOP for any 
station pair. Therefore, referring to the mathematical mo- 
del of Chapter II, an essentially equivalent problem is to 
fix L in (6a)and vary x in (6a)and (6b)over a range of values 
necessary to traverse all candidate LOP in the unambiguous 
lane given for the two stations being considered. Figure 6 
illustrates the geometry of the problem. This method is 
convenient for the following reasons. First, only one para- 
meter, x, need be varied in the estimation procedure; and 
second, for a fixed step size in an iterative solution, the 
total range over which x is varied (one unambiguous baseline 
lane width) is the same for any receiver position. A reason- 
able value to assign to L is the sum of the great circle dis- 
tances from the two stations to the center of the presupposed 
area. If х, is defined as the great circle distance from 
the reference station, henceforth referred to as station one, 
to the center of the presupposed area, the range of values 


over which x must be varied is х, + one-half the unambiguous 
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lune size. The fact that this range of values for x actually 
traverses an unambiguous lane is easily seen if one assumes 
two artificial stations separated by a baseline distance L. 
Estimation of x on the artificial baseline formed is equiva- 
lent to estimating x ot dde baseline in the real problem if 
the manner of estimation explained above is followed. Thus, 
the range of x which must be examined is one baseline unam- 
biguous lane width. 


To estimate x in the problem posed, we form the likeli- 


hobd function [26] as 


A(x) = f(Y/x) -JJ MOV а е СО Бе 0 (9) 


where V is the 2K-dimensional phase velocity column vector, 
0 is the K dimensional random phase vector and f(V) and f(9) 
are the joint probability density functions of the components 
of V and 0 which are assumed mutually independent.  f(Y/x) is 
the joint probability density of the components of the vector 
Y conditioned upon knowing x. Since our only interest in 
the likelihood function is to find the values of x which 
yield its maxima, we can equivalently maximize any monotonic 
function of it. Factors multiplying the likelihood function 
but not containing x, v, or 9, will be dropped without re- 
defining the likelihood function throughout the dissertation. 
Ambiguity resolution information can be conveniently in- 
corporated in the maximum likelihood estimator in the follow- 
ing way. Consider the ambiguity resolution information as 


a priori information about x in the form of a probability 
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density f(x). Equation (9) when multiplied by f(x) is pro- 
portional to the a posteriori density of x. The maximum a 
posteriori (MAP) estimate of x is then that x which maxi- 
mizes this equation and is the optimum estimate of x for 
the problem [26]. When f(x) is a uniform density function, 
as would be the case if exactly one unambiguous lane were 
given by an ey resolution technique with no informa- 
tion about where in the lane the correct LOP might be, the 
MAP estimator is equivalent to the maximum likelihood esti- 
mator. 

It can be shown [26] that for the density function of 


N(t) assumed for the mathematical model in Chapter II, 


Т 
2K 
£(Y/V,0,x)o 1 stamp зүсээ, будаа (10) 
qi 


о 
where T is the total data collection period. Substituting 
(4) in (10) and using the trigonometric identity for the 


cosine of.the sum of angles, (10) can be written as 


ЗІ 
о ыы ехр (5,8) L6; COS (6.791) 
+ S, sin САН (12) 


Ша (11) с. апа 5. are defined as 


~ 10(1-1) 4 Т, 
Cy 57, y, (t) cos[(w,+A,)tldt (12a) 


Ші 
IBID 
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И! 


2 


v 


and ШО = Ti 


M 
5, = 57 y,(t) sin[(o,*A,)t]dt — (12b) 


= 
: 10(1-1) 


where M is the number of T, sec. duration pulses contained 
in the data collection period. 

The quantities А.» N, and A; (1<i<2K) must be estimated 
from the collected data. The estimation equations for these 
quantities are derived in Appendix A and Appendix B. From 
this point in the two-station estimation derivation it will 
be assumed that the values of А. and N, have been estimated 
from the data and the estimates are available. The values 
of A. will be assumed to have been estimated and the data 


corrected for doppler shift utilizing this information. Fol- 


lowing these assumptions, we define 


үг” 2 


2 2 
1 і 十 Si ) (24,/N1) DS (13a) 


and 


tan ^ (S,/C,) 1<i<2K. (13b) 


те 
|| 


The likelihood function (9) may then be written by substitu- 


tion ot (13) into (M) into (9) as 


K 
A(x) = Ey J£ CO exp >, [Ricos(ó-u)tRi,,cos(ó!,,70,,,2]cos01 


1-1 
К 
ID ус А, үзіл (0, ada о, } 49 
1 = 


(14) 


where E, C? denotes expectation with respect to V. 


54 SN 





Since the 9, мене arsued tofbefindependent, uniformly 
distributed random variables on the interval (-T7, m] in 
Chapter II, we may use the relationship 

Ti 


L2 (15) 


(21) ^1 | exp[A cos Ө + В sin 0140 =1_[(A°+B°) 
-- ТІ 
some algebra and some trigonometric identities to write the 


likelihood function as 


K 
в 2 2 12 
ПОЈ Er, 2 мэн + ВО ВО 2В В „Со 28 } (16) 
where 
Завь LÀ » 1 en и 
Ui = principal value of TT ur ф 1194). (17) 


їп (15) апа (16) 1,0) is the zero order modified Bessel 
function of the first kind. At this point, no approximations 
have been made and the assumed random nature of the effec- 
tive phase velocities has not been utilized. Therefore, if 
the phase velocity column vector V is assumed to have deter- 
ministic onen с as was first thought to adequately model 
the problem, the maximum likelihood estimation equation for 
known phase velocities can be obtained. This estimation 
equation is derived in Appendix C. As was discussed in 


Chapter I, the phase velocities are not well enough known 





1тһе value of 0 which lies in the interval -п<б<т 18 


called the principal value of 0. 
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to use this model, so we continue with the derivation of the 
estimator assuming Gaussian phase velocities. 

A few minutes of data collection is expected to be ade- 
quate to insure high SNR in one Hz even for worst case noise 
[2]. The arguments of cha 1, function associated with the 
ith frequency are on the order of the sum of the SNR in one 
Hz of the signals received from the two stations at the ith 


frequency. This allows the replacement of the I. functions 


in (16) by their large argument approximation [27], 
І. (2) = a ийг ехр (2). (18) 


Using (18), the likelihood function (16) may then be written 


as 





K 
2 2 1/2 
З ГУ + Russ + 2R,R,,,Ccos Y] } 
A(x) - E • (19) 
CM 2 1/4 
ч + КЕНИ + Е | COS ф 1 


The numerator of (19) is a rapidly varying function of x, 
sharply peaked at all the local maxima of the function. The 
denominator is a more slowly varying function and for regions 
near a local maximum.can be considered to be nearly constant. 
The denominator will have different values from peak to peak 
of the likelihood function, but the magnitude of the change 
from peak to peak will be quite small compared to the change 


in the magnitude of the numerator for the same peaks. For 
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the peaks in the vicinity of the global maximum which we are 
seeking to find, the denominator in (19) can be considered 
to have relatively little effect on the likelihood function 


and will be ignored, yielding 


12) El 5 Гоа 2 т вв. cosy, 12 23. (20) 
2 үес? еі K+i і"кчі і е 
It is eonvenient to write 
v4 (x) = У; (>) + РҮ (21) 


where v; (x) is the expected value and p, some random pertur- 
bation of the ith phase velocity. The expected value of the 
phase velocities are all very near c [10], so we can approxi- 


mate (6a) and (6b) by 


x 


фто | и -2m 1<i<K (22a) 
5 у. (У) [1+Е, ] 
E 
жен ПН ү СЕБИНЕН 5 ; 
Pray 0. 27m, | x 1<i<K (22b) 


Ма СТО) ПЯ ча) 


where 


sj = p,/c - p,/v GO. (23) 


The random variables Ea» 1<i<2K will be taken as jointly 
Gaussian with zero means and as components of the 2K dimen- 
sional vector €. This follows directly from the assumption 


of the components of V being jointly Gaussian and the 


x i X 





FetfinitTon of €. The covardanee matrix oft с will be 





K (24) 


where K is the covariance matrix of the 2K dimensional 


vector V. Utilizing this result and the approximation 


T 





EET ч 1 - n for |а|<<1, we can now write 
E: BS 
22 = = <i< 
A = = №. 1<1<К (25a) 
and 
01 (L-x) 
а тои eee E е 
where 
х 
Ф; = «| = 1<1<К (26а) 
UD 
о 
and 
L 
E и ed Щщ | 
2277 ш, = К! (26Ъ) 
xo KH y 


Recalling the definition of МА from (17), for high SNR 
Vi will be small near the local maxima of the likelihood 


function and we write 
2 
i 
2 e 


р 





cos E ME (27) 
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This allows the likelihood function to be written as 


2 
Карт 1/2 
LA C 
) 


K 
ЕЕ exp { 2 ФЕ + Е )[1- 
€ bd K+i (R,*R 


. (28) 
Kt i 
Again, using the fact that Ui and thus у. is small in the 

vicinity of the local maxima of the likelihood function and 


the approximation 


пря 


(1 - п) W= n2, (29) 


an.equivalent likelihood function can be written as 


K 
u 2 
A(x) = Е, ехр [Е 0,8, 1 (30) 
ї=1 
where 
R.R : 
Q, T K+i ‚ | (31) 
2 (RR Li) 


Using vector notation, the likelihood function can be 


written as 
Т T 
A(x) = E. exp [e "Ae + Be + Р]. (32) 


The 2K-by-2K matrix A is defined by jth row and nth column 


entries 

0) 2 “50 
a, = — l«j-2nc«kK (33a) 
jn 2 = = 

с 
Ті ex Q, K 
EON жен |а 5, и 

m 5 K+1 < j n < 2K (33b) 
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Ben. у = 
а ја 2 1 < j =|(n+K)| < 2к (33c) 


Зам = 0 elsewhere. (33d) 


The 2K-dimensional column vector B is defined by ith entry 


bi = - 2 910; х П к (34a) 
bi - = E КЕГІ < i « 2K (34b) 
where 
Ч, = principal value of el: (35) 
The scalar D is defined as 
Ч 2 
D = d и, . (36) 


The expected value in (32) under the conditions given is 
shown in Appendix D to yield 


DDr 2AK ) 5B] 
А(х) = 8—2 (37) 


|І - 2AK 122 
E 


* 


denotes the 








where I is the 2K by 2K identity matrix and 
determinant of a matrix. The value of x which maximizes 


the likelihood equation can be determined by allowing x in 
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(37) to range over the necessary values as shown in Figure 6. 
Other less likely LOP can be determined and a measure of the 
probability of their being correct can be found by dividing 
КО У атеке Еос the Bess) likely LOP by the value’ of A for 
the most likely LOP. If the LOP being estimated is considered 
to be a random variable, uniformly distributed across one un- 
ambiguous lane, this ratio is a ratio of the a posteriori 
-probability densities of the two LOP. 

Figure 7 is a drawing of the melone of a typical likeli- 
Mood function as it varies with x. Within the "envelope" 
the function will be a more rapidly varying function of x. 
The number of local maxima in an unambiguous lane, determined 
by the values and number of Omega frequencies, will in 
general be many more than shown. The value of х,» the dis- 
tance from the first station to the center of some presup- 
pased area, serves only asia means to establish the fixed 


value of L to use in the estimation. In Figure 7, х the 


ml’ 
value of x which yields the maximum value of Á cannot be 
used as a valid measure of the distance of the receiver from 


the reference station, but serves only to identify the re- 


sulting LOP estimate. 


В.  TWO-STATION ESTIMATION FOR FOUR FREQUENCIES 

The general estimation equation for K frequencies re- 
ceived from two stations is given by (37). This equation 
can be applied to the four-frequency case with the size of 


the matrices A and K being 8-by-8 and the column vector B 
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being 8 dimensional. For the four-frequency, two-station 
problem then, the likelihood function is given by 


exp[D*(1/2)B K, (I-2AK ) ^ 8] 
Suy GALL жиса ээс (38) 


B2 
pg 


where К is the 8-by-8 covariance matrix of the normalized 
effective phase velocity vector €. The 8-by-8 matrix A is 


given by 
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(6€) 

,(к-1)70, 7 0 

0 г (жопу, бо 

0 0 

0 0 
(x-1)x"o, "n- 0 

0 (x-1)x°D,, m- 

0 0 

0 0 


г 


Ск-т)ж7о m- 


0 


(хет) x*b, *n- 


(x-1)x, b, "те 0 


2 
т. IT С 
X 0- ——— J 
0 (х-т)х 0, I 
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where Q,, I u emo given by (31). The column vector ® 


in (38) is given by 


xu 
205 ost 
E E 
-0,М,ш,Х 
B = = 898? (40) 


919% (0-х) 
0,0,0, (1-х) 
0503 (2-х) 


Q,W,0, (L-x) 


where Wi, І < і < 4, іс рімет by (55). — Ihe scalar D is 


given by 


站 二 NE 3 te (41) 


C. ESTIMATION ERROR TYPES 
The cyclic nature of the likelihood function yields LOP 
estimates which fall into three categories: 
1. Estimates which fall within nominally 1 or 2 nmi or 
less of the correct LOP, 
2. Estimates which are nominally one half wavelength 
of the fundamental four frequencies from the correct 
LOP, referred to as "minor lane" errors and 
3. Estimates which are farther from the correct LOP than 
the previous two types, which we refer to as "major 


line" errors 
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The latter two types of estimates are considered as erroneous 
estimates for the purposes of this thesis. An important con- 
sideration to note is that the probability of making a posi- 
tion error when using data from several stations will very 
likely be less than the probability of making a lane error 

in the two station case for the following reasons. The trans- 
missions of at least three Omega stations will be available 
(this number being required for a position fix) and it is 
very likely that one or more additional stations will be 
usable given the great distances of VLF propagation and the 
planned worldwide location of eight transmitters in the Omega 
System. This redundancy will reduce the probability of large 
position estimation errors. The multiple station problem is 


considered in Chapter VII. 
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V. FOUR-FREQUENCY SIMULATION 


In the interest of achieving unambiguous lane widths 
which can be identified by an ambiguity resolution technique, 
additional frequencies may be added to the Omega format. 
Assuming that the present three Omega frequencies are re- 
tained, adequate widening for ambiguity resolution by time- 
of-arrival techniques [21] is expected to be achieved by 
Madi tion of a fourth frequency. A fifth frequency and/or 
frequency modulation of one of the signals by a signal with 
wavelength on the order of several thousand miles would 
yield better estimator performance and perhaps inherent 
ambiguity resolution as well. This chapter considers the 
performance of the four-frequency Omega maximum likelihood 
estimator for several candidate fourth frequencies which 


can be added to the Omega format with minimum modification. 


A. FREQUENCY CHOICE 
The criteria used in selection of the fourth frequencies 
were: 
l. су ш: the present three Omega frequencies 
20.2 kHz, 11 1/3 kHz, and 13.6 kHz), 
2. Separation of the fourth frequency from any of the 
existing three Omega frequencies or other known 
existing or planned VLF transmissions by at least 


250 kHz to avoid interference in present receivers, 
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3. A resultant unambiguous lane size greater than the 
present three frequency baseline size of nominally 
72 nmi and 
4. Selection of a fourth frequency which can be 
generated using the present Omega signal format 
generator frequency of 816 kHz. 
In addition to these criteria, preliminary four-frequency 
simulations were carried out for SNR which are to be 
expected for reasonably long transmission paths and assun- 
ing that exact phase velocity predictions were available. 
These simulations, discussed in Appendix C, indicate that 
ee unambiguous lane widths must be restricted 
to about 1000 nmi or less for estimation performance "accept- 
able" for SAR application. The frequencies which met the 
four criteria above, the integer divisor of the fundamental 
station frequency (816 kHz) and the resultant nominal base- 
line unambiguous lane width are presented in Table III. 
Figure 8 shows the Omega frequency spectrum with the present 


three frequencies and the potential fourth frequencies. 


В. SIMULATION DATA AND PROCESSOR 

Monte Carlo simulation of the four-frequency Omega prob- 
lem was carried out utilizing the maximum likelihood esti- 
mator given by (37) in Chapter IV with K = 4, The parameters 
which were varied were the choice of fourth frequency, the 
SNR in one Hz and the uncertainties of the phase velocities. 


It is shown in Appendix A that C, and 5.» when conditioned 


i 
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TABLE III 


CANDIDATE FOURTH FREQUENCIES 


Nominal Unambig. 


Frequency > Divisor (n,) Lane (1/2 A) 
10.4615 kHz 78 928 nmi 
' 10.8800 kHz 75 357 nmi 
11.6571 kHz 70 500 nmi 
69 
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upon knowing M and m and also under the assumption that 
N(t) is a vector whose components are sample functions of 
a white, Gaussian process, are Gaussian random variables. 


The conditional means and variances of Ci and 5. мһеп 


е 


doppler shift is assumed to be negligible are given by 


A ML. 


Е (с. |ф. „0. ) En cos QUUM T 9.) : (42а) 
А МТ. 
Е (5. |. „6. 1 EN (42b) 
and 
N MT, 
VAR |С,19,»9,1 = VAR [5;1Ф;› И (42 с) 


For the simulation, A; and 9. were assumed equal to zero and 
А. and N, were determined from the SNR in one Hz, assumed 
the same for all i. The product MT, was assumed equal to 

18 for all i, this choice being influenced by the GRAN anti- 
cipated 180 second data retransmission time, equivalent to 
18 nominally one second pulses of Omega data for each fre- 


quency from each station. With these assumptions, the means 


and variances of с. апа 5. conditioned upon m and 9. become 


Е [С |ф,, Bl DA E (43a) 


E 15.19,» то (43b) 
and 


ON. 
i 


2 e 





VAR ІС,|%.. PNE УАК15, 19,» Ја = 


(43c) 


The eight phase velocities used to generate the data were 


generated by a correlated Gaussian random number generating 
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computer program using the mean vector with all components 


equal to c and the covariance matrix 


колесо : (44) 


In (44) p is a symmetric 4 by 4 correlation matrix referring 
to the interfrequency correlation of transmissions from each 
station and o represents a 4 by 4 zero matrix indicating un- 
correlated interstation phase velocities. In (44) с? is the 
variance of all Е; (equal to the variance of у; / с) which is 
varied as a parameter in the simulation. The correlation 
matrix was chosen to reflect uncorrelated interstation trans- 
missions because in the general case very low correlation is 
expected [19]. The inter-frequency correlation matrix p was 
formed by assuming the correlation of 10.2 kHz and 13.6 kHz 
frequencies from one station to be .7, which seems to be a 
reasonably good choice [16, 18, and 19] and linear interpola- 


tion using the equation 


Pij = 1 - 34309 | iT fp) Ca) 


where P1 4 isi thei ith row, jth column entry in the matrix® 6 
and Е, is the ith frequency. These assumptions resulted in 


the symmetric matrix 
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1.00 „90 270 014 
. 90 1.00 „80 Da 
MD „70 „80 1.00 Ди 222 
ВА P42 P43 RU 


where the first three rows/columns of the matrix correspond 
to the frequencies 10.2 kHz, 11 1/3 kHz, and 13,6 kHz res- 
pectively. The values used for E Pi? 1 = 1, 2, 3 were 
chosen from Table IV depending upon the fourth frequency be- 
ing used. 

The values for L and x used in the creation of the simu- 
lation data were 8000 nmi and 3000 nmi respectively. In the 
final Omega system implementation, these values should be 
typical of values to be expected. With all the above assump- 
tions, the simulation data was generated by computer using 
random number generating programs for the random variable 
sample values. 

The processing of the data with the four-frequency maxi- 


mum likelihood estimator, was accomplished using c as the 


mean value for all phase velocities and the covariance matrix 





K (47) 


where Ку is given by (44) and the value of the parameter or 
Using c for the mean value of phase velocity for processing 
is equivalent to having a reasonably good predicted mean 
value of phase velocity since the data were created using 


this mean value. 
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TABLE IV 


FOURTH FREQUENCY CORRELATION COEFFICIENTS 


Frequency (£,) Ч 0:4 
10.4615 КН2 1 ‚98 

‚92 

272 

10.8800 kHz 1 ‚94 
‚96 

3 2706 

he Ok HZ 1 "87. 

. 97 

3 ‚83 
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Figures 9, 10, and 11 illustrate the rekults of the simu- 
lation for each of the three choices of fourth frequency. 
Each point plotted consists of a number of "runs" which cor- 
respond to a three minute data collection period. The frac- 
tion of the total number of runs for each point which yielded 
lane errors (either major or minor, as defined in Chapter IV) 
ame plotted as a function of SNR in one liz and 0, the standard 
deviation of the phase velocities. Table V gives the number 
of runs which were utilized in the generation of each point 
and a breakdown of the number of the lane errors which were 
major and minor. 

The ERE, of runs which were used in generating the 
points in Figures 9, 10, and 11 are quite small for simula- 
Bion, particularly considering the small number of resultant 
Vane errors. The reason for this is the computer time which 
is required, 15 to 40 seconds on the IBM 360 per run depend- 
ing on the lane size. However, the reasonably consistent 
increase in lane errors with increasing lane width and phase 
velocity standard deviation influences the decision that 
the general 211: аге уа114, 

The results of this simulation prompt the conclusion 
that with adequate knowledge of the phase velocity variance 
and correlation and with reasonable signal-to-noise ratios, 
the unambiguous lane size can be increased beyond the present 
nominal 72 nmi obtained with three frequency Omega. The 
controlling factors in the selection of a fourth frequency 


for Omega to be used for SAR purposes are (1) the capability 
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SIMULATION RESULTS 


TABLE V 


Lane Errors 





Frequency (£) Std. Dev. (o) SNR No. Runs (major)(minor) 
10.4615 kHz 2 x 1074 10 dB 100 0 
" 15 dB 100 
" 20 dB 190 3 2 
вв 1072 10 dB 100 ey, 16 
" 20 dB 100 17 
10.8800 kHz тос" 10 ав 100 
" 20 ав 100 
же" 10 dB 100 
E 15 dB 100 
" 20 dB 200 
3 x n2 10 dB 100 18 
" 17 dB 100 20 
" 20 dB 200 0 26 
11.6571 kHz l x 1074 10 ав 100 
" 20 dB 100 
2 х 107“ 10 ав 100 
" 15 dB 100 
" 20 dB 200 
3x 1074 10 ав 100 3 17 
" 20 dB 100 2 10 


таз: 





of ambiguity resolution techniques to resolve the lane am- 
biguity with restricted amounts of Omega data and (2) the 
accuracy with which the phase velocities of the Omega fre- 
quencies can be predicted. 

One of the results of the Monte Carlo simulation as 
shown by Figures 9, 10, and 11 is that as unambiguous lane 
size is increased, the probability of lane error increases. 
However, the estimator provides several estimates along with 
their relative probabilities of being correct. Thus, even 
Em a lane error occurs, the second or third most likely 
estimate is nearly always correct. This fact, in conjunc- 
Биоп With) thesweduced probability of position error expected 
when redundant information from additional stations is avail- 
able, makes a .1 or .2 probability of lane error more 
acceptable as a trade-off for the larger unambiguous lane 


size. 
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УГ. А FTOUR=PREQUENCY OMEGA EXPERIMENT 


All of the ambiguity resolution techniques discussed in 
Chapter III, which appear practical to implement, require an 
unambiguous lane width in excess of the nominal 72 nmi which 
the present three Omega frequencies yield. The performance 
which can be expected from a maximum likelihood estimator if 
used with four-frequency Omega for several easily implemented 
eR frequencies was presented in Chapter V. Naval Air 
Test Center, Patuxent River, Maryland conducted a four-fre- 
quency experiment utilizing one of these frequencies (10.88 
kHz) in the fall of 1973. The maximum likelihood estimator 
derived in Chapter IV was utilized to process this data, the 


results of which are reported in this chapter. 


A. DATA COLLECTION 

The four-frequency Omega data was collected by Texas In- 
struments Corporation at Dallas, Texas [29]. The receiver 
location was 1197.53 nmi from the Forestport, New York experi- 
mental Omega station and 2349.62 nmi from the Trinidad Omega 
station. The data was collected over the two five-day periods 
17-21 October 1973 and 14-18 November 1973. Data was en 
for three-minute intervals four times per hour beginning on 
the hour and half hour and seven minutes past the hour and 
half hour for the entire periods. Each of these three-minute 


collection periods, herein referred to as "runs," consisted 
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of 18 nominally one-second pulses of signal for each frequen- 
су from each station. Оп the last day of the first data col- 
lection period, the Forestport station reduced transmitted 
power to very low levels. This portion of the data is not 
included in the estimation results as it is considered not to 
be representative of actual conditions to be expected. 

The receiver used for the data collection was the Omega 
Position Locating Equipment (OPLE) ground station modified to 
receive four frequencies [29]. A block diagram of the receiv- 


ing equipment is shown in Figure 12. 


В. MAXIMUM LIKELIHOOD ESTIMATION 

The maximum likelihood estimator for high SNR Omega was 
derived in Chapter IV. Quite high SNR were experienced as 
shown in Appendix B. As described in Chapter II, the effec- 
tive phase velocities are modeled as being jointly Gaussian. 
The eight element column vector V is defined as the effective 


phase velocity vector 


T 
У = [v,; Vos e e vgl : (48) 


The first four elements of V correspond to the four signals 
ЕН fF nequencies 10.2 kHz, ІІ 1/3 kHz, 13:6 КНЕ, апа 10.88 kHz, 
in that order, from Forestport. The last four elements of V 
correspond to the same frequency signals from Trinidad. The 
mean values of the Gaussian vector V which were used in the 
estimation were those provided by Professor J. A. Pierce of 


Harvard in [16] and in private correspondence. These values 


82 





е долтеоещ цот209ж11029 EJLA “721 елп та 


0 


3 0 UTS 


| 


|) 
f^ 
И 


)) 
| 


= 
f\\ 


(0 500 


я 
[т] 
Mes СУ 
Е 
us quo Е < 
Я4У1 В Ы = = 
TVII9IQ = e -% Ч 
© . 2 
БЕ > ---- 5 је 
= : 1 = 
ЁО » - = — 
S 5 
4 Q2 





83 





are given as ratios of the phase velocity and the velocity 

of light in vacuum (v/c) in Table II. Because of the signi- 
ficant diurnal effects of phase velocity, different values 
for v/c are required for day (all lighted transmission path) 
and night (all dark transmission path). For transition per- 
lods, when a path was not all lighted or all dark, weighted 
averages of the day and night path phase velocity values from 
Table II based upon the fraction of the path which was light 
was used. In general, this would not be justified, but since 
the values of phase velocity for day and night are so close 
and because of the phase velocity uncertainties, it is con- 
sidered acceptable to use this technique. 

Since the phase velocities are considered to be jointly 
Gaussian, a covariance matrix must be specified in order to 
do the estimation. The inter-station correlation of phase 
velocities was TE еен to be zero and the inter-frequency 
correlation was determined by assuming the phase velocities 
of the 13.6 kHz signal and the 10.2 kHz signal to have a 0,7 
correlation coefficient. These assumptions, taken from [16, 
18, апа 19], allow the matrix О to be determined as in Chap- 


ter V, Equation (46) 


1.00 „90 «70 „94 
ЕЛІ 2.00 .80 Eoo 
. 70 . 80 1.00 ЕЛО 
. 94 . 96 „76 1.00 
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The matrix p is used, as in Chapter V, to form 


КА с (50) 


where c is the speed of light in vacuum, О 15 the standard 
deviation of v/c which was assumed to be the same for all 
frequencies, and о is the 4-by-4 zero matrix. The final 
parameter which must be specified for the estimation pro- 
cedure is О. As discussed throughout this thesis, the value 


of v/c is expected to be predictable to a few parts in Ножі 


2 


Values of o from 10 ° to 5 x To were used in the estima- 
tion presented here. It is expected that future analysis 
of Omega data for these frequencies with the idea of con- 
Sidering predicted phase velocities will yield better values 
for р апа О. It is not unlikely that different standard 
deviations of predicted phase velocities will result for 
different frequencies. If this does occur, Ку їп (50) сап 
easily be modified. 

with the assumptions delineated above, the estimation 
can be performed. As discussed in Chapter IV, the estima- 
tion procedure consists of assuming a priori knowledge of 
an unambiguous, nominally 360 nmi, lane which contains the 
correct receiver LOP. This information must be provided by 
an ambiguity resolution technique in a final system imple- 


mentation. Then, the sum of the great circle distances from 


the two stations to some point in the region of the receiver 





is fixed equal to L. Since two stations can улеа only an 
LOP, this sum of distances can be quite approximate, as will 
be the case, since this sum must actually be obtained from a 
coarse ambiguity resolution technique. Finally, x is de- 
fined as the distance from a reference station (Forestport 
was chosen for this processing) to a potential receiver site. 
Evaluation of the likelihood function as x is varied over 
the known unambiguous lane provides a function which has its 
e mimi value £r thatax which is thel maximum likelihood 
estimate., This estimate of x (x 1? is equivalent to estima- 
tion of an LOP as can be seen from Figure 6. The path over 
which the evaluation is done is shown by the dashed line, a 
convenient path for the following reasons. First, only one 
parameter, x, need be varied in the estimation procedure; 
and second, for a fixed step size in the evaluation, the 
total range over which x is varied (one unambiguous baseline 
lane width) is the same for any receiver position. This 
method of estimation necessarily yields all of the local 
maxima of the likelihood function. Thus a selected number 
of estimates of x which are less likely to be the correct 
estimate than x j сап easily be retained. As discussed in 
Chapter IV, the ratio of the likelihood function evaluated 
at these alternate ар а values of x divided by the like- 
lihood function evaluated at x 1 is a measure of the relative 
probability of the respective values of x being correct. 
Figure 13 shows a block diagram of the estimation pro- 


cedure. The collected Omega data is provided to a processor, 
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МЕ со ред за дева! 1 in Chapter IV, which evaluatesja likeli- 
hood function over a prescribed range of values for x. Para- 
meters which are provided to the processor are the mean 
values and covariance matrix for the phase velocities. 

The IBM 360-67 at the Naval Postgraduate School was 
utilized for the processing of the Omega four-frequency datas 
A Fortran program (Appendix E), requiring 34K BYTES of core, 
was used to iteratively evaluate the likelihood function for 
x at 0.5 nmi increments from x = 855 nmi to x = 1540 nmi. 

The one nmi regions centered at the values of x which yielded 
the aioe largest values for the likelihood function were 
then iteratively searched with 0.05 nmi increments to obtain 
a fine estimate. The values of the likelihood function for 
the second and third most likely values of x were then 
divided by the value of the likelihood function tox х 1? 
yielding a probability ratio discussed earlier. The time 
required to process one run (a data set corresponding to 
three minutes of collected data) was on the order of 15 
seconds of computer time. The range of x for which the pros 
cessing was done was nearly twice the nominal 360 nmi which 
would have been required had an ambiguity resolution tech- 


nique provided the unambiguous lane over which x was to be 





Any number of values may be retained, however, three 
were found to be adequate to always contain the correct 
distance. | 


> 
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varied. This was considered the only fair way to process 
the data since the correct value of x (1197.53 nmi) was 
known a priori and was considered as the center of the un- 
ambiguous lane. In a finally implemented system, the a 
priori knowledge of x will consist of some probability den- 


sity of x, across one or more unambiguous lane widths. 


ESTIMATION USING PREDICTED FIRST MODE PHASE VELOCITIES 

As an intermediate step, the optimum estimator must esti- 
mate the amplitudes, noise density heights and doppler fre- 
quencies of all eight signals. Appendices A and B give the 
derivation of these estimators. The doppler shift of the 
signals, which would result from physical motion or frequency 
drifts in the receiver, were found to be negligible since 
there was no SARCOM and thus no motion and were assumed to 
be zero for this experiment. In a finally implemented systen, 
doppler Er: must be estimated and removed. The signal 
amplitudes and noise density heights which were estimated as 
an integral part of the processing were used to form the sig- 
nal power-to-noise density ratios in dB-Hz. The results 
of this estimation are presented in Appendix B. 

The ec of the LOP estimation utilizing the optimum 
estimator are presented in Table VI. The mean values of 
phase velocity shown in Table II, and the values for Ку given 
by (49) and (50) were used in the processing. Several values 
of 0, the standard deviation of the phase velocity to velocity 


of light (v/c) ratios, were used, but only O = 2 x ЛОЛА о = 
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3x ови. апа с = 4 х о are shown. The results of using 
other values of O were significantly worse than for those 
shown. The minor lane errors which are referred to in Table 
VI are defined as runs for which one of the peaks of the 
likelihood function adjacent to the correct peak is the lar- 
gest value of the function. This results in an LOP estima- 
tion error which is nominally 6-7 nmi from the correct LOP 
on the baseline. Major lane errors would be defined as lane 
errors greater than a minor lane error. However, in this 
experiment using the values of О stated above, major lane 


"per- 


errors did not occur. In Table VI, the column labeled 
iod" okei the numbers 1 ог 2 апа the letters D, N, or T. 
These refer to data collection periods 1 252 and Day A Night; 
or Transition condition paths, respectively. 

The means and standard deviations of the LOP estimation 
error are in baseline measure., The hyperbolas which are 
formed by the isophase lines from a station pair in a phase 
comparison system diverge off the baseline (great circle pass- 
ing through the two stations). The divergence occurs directly 


' where 8 is the angle 


as csc (0/2), ene “spreading factor,’ 
formed by the two great circles passing through the receiver 
location and the stations as discussed in Chapter I.  There- 
fore, phase differences at points off the baseline correspond 
to greater distances than on the baseline. The spreading 
factor for the Dallas receiver location for this experiment 


was 1.88, thus the distances in Table VI must be multiplied 


by 1.88 to obtain the true measure for Dallas. The practice 
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of using baseline measure was selected so that the results 
сан вейшате readily аррітей со rhezgecneral case where a 
spreading factor of from 1 to perhaps 2 might result. 

The values of o which yielded the fewest errors were 
different for the two das collection periods. The value 
of o = 2 x 10 resulted in the fewest minor lane errors for 
the first period and 0 = 4 x 1074 resulted in the fewest for 
the second period. Increasing the value of oO is equivalent 
to desensitizing the estimator to the mean values of phase 
velocity being used. It would thus appear that for the 
data obtained in the second period, the mean values or pre- 
dicted ic of phase velocity which were used were not as 
close to the actual phase velocities which were encountered. 

It is believed, however, that the data was significantly 
disturbed by higher order modes of propagation for the night 
path from the Forestport station. Therefore, use of the 
first mode phase velocity values given in Table II would not 
be expected to yield error-free estimation results for 
night paths. This conclusion was reached after the follow- 
ing considerations. The mean values of the phase differences, 
Ab, for all four frequencies, f, for all ТЭЭХ (day) trans- 
mission paths were utilized to estimate v/c from the data 
by using 

f. (L-2x) 


(v/e); = —— 1<1<4 (51) 
e(n; + Аф. /2т) 





to obtain Figure ila. ico) eee ice Chie, (Known) sumfof the 
distances from the stations to the receiver, x is the (known) 
distance from the reference station to the receiver, c is 
the velocity of light in vacuum, and ny is an integer equal 
to the number of whole cycles of phase difference at the re- 
ceiver. Equation (51) assumes that the phase velocities for 
the same frequency signals were the same for both paths under 
day conditions. This is considered to be reasonable since 
the paths were both predominantly westerly and first mode 
propagation during day conditions apparently dominated. This 
latter assumption seems to be borne out by the fact that no 
day TU lane errors were experienced. The apparent 
parallel shift of the curves in Figure 14 from the predicted 
first mode values could be explained by timing errors be- 
tween the two stations. Several 0.25 yus and 0.5 ys correc- 
tions in timing were inserted by Forestport during the data 
collection periods. This parallel shift in phase velocities 
will not cause lane errors to result, but will be seen as an 
estimation error or a bias in the estimation error mean value 
for this 2 755711) 

The results of Figure 14 can now be used to obtain phase 
velocity values for dark paths (night) for all four frequen- 


cies from each station separately to determine if the phase 


EE predicted phase velocity values shovn in Figures 14, 
i, 16, and 17 are taken from (11, 12, and 16] and private 
communication with Professor J. A. Pierce of Harvard Univer- 
sity. The lines drawn between the plotted points serve only 
as aids for identification of points and should not be con- 
strued as portions of continuous curves. 
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velocities from one station appear to be causing the lane 
errors. This was done by using the change in the mean values 
of the phase differences from day to night for each frequency 


from each station separately in the equation 


ДЕ е): 
(wie), ya — 16444 . (52) 


2 c(i N - PT RUE, + 2nf,d 


In (52) d is the distance from the station under considera- 
tion to the receiver, (v/c)i p are values taken from Figure 
3 


14, and (4. ) is the difference of the mean values 
1,N 


Фі 
of the night phase and the day phase for the ith frequency, 
treating the stations separately. The results of the night 
phase velocities obtained in this way are presented in Fig- 
ures 15 and 16 for the first and second data collection 
periods a ee The phase velocities from the Forest- 
port station are, by this analysis, apparently the cause of 
the lane errors which were experienced for the night condi- 
tions. The relatively parallel shift of the Trinidad phase 
velocities from the predicted first mode values indicates 
that the Trinidad data is probably dominated by first mode 
Propagation. However, the erratic shifting of the Forest- 
port phase velocity values from the first mode predictions 
indicates very probably that higher order propagation mode 
interference was experienced during the night conditions for 


both data collection periods. 


Ibid. 
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One further processing of the data was performed to 
demonstrate the consistency of the estimator. Phase velocity 
values for each frequency for lighted transmission paths and 
for each frequency for dark transmission paths were deter- 
mined for each data set б (5): This method determines 
effective phase velocities, assumed the same from each sta- 
tion for the same frequency and path conditions. The phase 
velocities, normalized by c, thus determined are shown in 
Figure 14 for lighted paths and Figure 17 for dark оо ил 
The median values of phase velocities of the two collection 
periods were then taken from the curves and the entire data 
set was processed using this one set of phase velocity 
values. 

The results of this processing were: three minor lane 
errors for the second period with night conditions and one 
minor lane error for the second period with transition condi- 
tions. All other data were lane error-free. This demon- 
strates that if reasonably accurate effective phase velocity 
values are known, extremely good estimation results can be 
expected from the maximum likelihood estimator. It further 
ELpeOnstradte's Ена ар least ior vhegexperiment conducted, the 

| 
mean values of the effective phase velocities did not change 
significantly in more than one month. A single phase veloc- 
ity value for each frequency for night and likewise for day 


conditions yielded excellent results for data collection 


periods separated by a month. 





Ibid. 
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As stated before, higher order mode interference is the 
suspected cause of the lane errors experienced in processing 
the data utilizing the phase velocity values from Table II. 
This and the demonstrated relative insensitivity of the esti- 
mator to phase velocity values, along with the lane error- 
free day path condition results, indicate that very good 
estimation results can probably be expected from using single 
phase velocity values for a given receiver location. The 
Бепейсіуісу of these phase velocity values to receiver loca- 


tion remains to be determined. 


D. ESTIMATION USING PROPAGATION CORRECTION TABLES 

The Defense Mapping Agency publishes Omega Propagation 
Correction tables regularly for the active Omega stations in 
pairs for the two frequencies 10.2 kHz and 13.6 kHz. These 
corrections, when added to the observed phase difference of 
the appropriate frequency from a station pair, allow a map- 
ping to be performed from the observed phase difference to 
an Omega navigational chart. These charts are plotted for 
an assumed value of (v/c) » 1.00261 for all frequencies. 
Many factors in addition to diurnal effects which affect 
phase, or equivalently the effective phase velocity, of the 
signals is accounted for by the computer program which gen- 
erates these tables. In order to determine whether the pre- 
dicted corrections would improve the results of the experiment, 
tables for the appropriate location, time, and station pair 
were obtained from Defense Mapping Agency. The corrections 


used for the 10.88 kHz and 11 1/3 kHz signals were obtained 
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by linear interpolation from the published corrections. 
Curves showing 24-hour plots of these correction factors 
are shown in Figures 18 and 19 for the second half of Octo- 
ber 1973 and in Figures 20 and 21 for the second half of 
November 1973. | 

The results of the estimation when the data were cor- 
rected using the factors found as described above and using 
(v/c) = 1.00261, for all frequencies is presented in Table 
VII. The number of lane errors was decreased by utilizing 
the propagation correction factors; however, quite a large 
number of lane errors still occurred at night. It is be- 


lieved that this helps to substantiate the higher order mode 


interference theory of the last section. 


Е. EXPERIMENT ESTIMATION SUMMARY 

Because of the suspected higher order mode interference 
at night during this experiment, it cannot be considered 
conclusive that use of the Table II phase velocity values 
are inadequate for LOP estimation in a final SAR system. 
Further experiments at points an adequate distance from the 
stations to preclude higher order mode interference may show 
that propagation correction factors are helpful in reducing 
lane errors. If this results, the corrections can quite 
easily be implemented at the SARCEN and are not considered 
as a major hindrance in implementing a successful system. 

One of the most important beneficial factors in using 
the maximum likelihood estimation procedure is the avail- 


ability of multiple estimates along with their associated 
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relative probabilities of being correct which are obtained. 
When lane errors were made using the phase velocity values 
Brome аъ ер for the first period, the correct LOP was 
always the second most likely. For the second period, the 
correct LOP was зэм ER likely for all but eight runs 
where lane errors occurred. For seven of the eight runs the 
correct LOP was third most likely and for one run it was 
fourth most likely. When propagation correction tables were 
used, the correct LOP was second most likely for all but one 
run where lane errors resulted. For this run the correct 
LOP was third most likely. 

For thie estimation performed with the Table II phase 
velocity values and © = 3 n оныг the probability ratios dis- 
cussed in Chapter IV were examined. As explained in Chapter 
IV, the ratios of the likelihood function evaluated for two 
values of x is equivalent to a ratio of the probabilities of 
these x being the correct x for the receiver location.  Fig- 
ures 22-25 give selected probability ratios for the first 
data collection period and Figures 26-31 give selected pro- 
bability ratios for the second period. The symbol 9 denotes 
thel probability ratio of the second most likely x when the 
correct x was the most likely. The symbol À denotes the 
probability ratio of the second most likely x when the cor- 
rect x was second most likely. The symbol [+] denotes the 
probability ratio of the correct x when it was third or 
fourth most likely, the latter case occurring for only one 


Тип ber en tre ага сен. 
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From these results, it appears that a threshold could be 
established such that when a probability ratio of an x was 
above a threshold, it would be considered as likely and SAR 
forces could be sent to the “үнд It is not suggested that 
the probability ratio values be treated as actually having 
meaningful numerical value. It appears that the maximum 
likelihood estimator places quite a lot more confidence in 
the data for this experiment than the results of the experi- 
ment would suggest is reasonable. Whether these values will 
T more meaningful as numerical values when receiver sites 
where higher order mode propagation is not a factor must be 
determined by future experiment. 

The information of "ranked" position estimates which can 
be provided is considered to be of great importance for a 
SAR system. If the most likely position estimate does not 
yield a recovery, rather than spending hours or even days in 
an ever-widening search, the SAR forces can be vectored to 


the next most likely position. 
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VII. MULTIPLE STATION POSITION ESTIMATION 


А. AN OPTIMUM POSITION ESTIMATOR 

The mathematical model to be used for the multiple sta- 
tion problem will be the same as that for the two station 
problem presented in Chapter II except that the vectors will 
be of length PK instead of 2K where P is the number of Omega 
stations from which data is received. Arbitrarily ordering 
the stations as one through P, the first K components of the 


received column vector 
YCE) sa Sty taN Ct) (53) 


will be from station one, the second K from station two, etc. 
Here again, K is the number of different frequencies common 
to all Omega stations for use in hyperbolic navigation. In 
(53), S(t) and N(t) will be PK dimensional column vectors as 
well. The components of S(t) will represent the signal por- 
tion of the received vector Y(t) and the components of N(t) 
will e Eau d functions from white, Gaussian noise 
processes, the ith component n; (t) having two-sided noise 
density height N./2.. The ith signal component can be written 
as 


s,(t) E А. cos [Co + A.)t - $1 + 9] (54) 
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as in Chapter II. The same arguments as followed in Chapter 


Il allow the following observations 


i КҮТ I II (РК (55) 


and 


9, = 0 = 0 22... 


і Кра 2K+i БА? 


9 (р-1)к+і 
where 0, CITA атс inde pendent) random ас тар ебе, uni- 
formly distributed on (-1, 1]. 

It will not be possible to follow the convention in 
Chapter II for defining the propagation phase delays.  There- 


fore, the following definition will be made 


“ын 


2 к аи . 
2701 = г vies, 09) un m (1<i<K) (57) 


a [o<p<(P-1)] 


where Vk is the kth phase velocity and X jis the- gneat circle 
distance from the nth Omega station to Eu receiver location 
being estimated. і 
Similar to the presupposed unambiguous lane in Chapter 
IV, a presupposed area on the globe is required to begin the 
multiple station estimation process. This must be obtained 
by some ambiguity resolution technique using the data re- 
ceived from P stations. The estimation problem now consists 


of evaluating a likelihood function for all points on the 


presupposed area and determining the column vector x which 


120 Ж 





maximizes the function on this area. The likelihood func- 
tion will appear as a rubber sheet spread over the pre- 
supposed area with many spikes under it causing many "hills" 
in the sheet. The value of the column vector x which cor- 
responds to the "highest hill" will be the maximum likeli- 
hood estimate of x and the components of this x will 
uniquely determine the great circle distances from the Omega 
stations to the maximum likelihood estimate of the Omega 
receiver position. 

Following the notation and arguments of Chapter IV, the 
likelihood function for the multiple station problem can be 


written as 


К Р-1 
A(X) = ЕС ERBE (кәні <6 ei 
і-1 р-0 
К Р-1 
cos 9, + 2 i EC sin (Ф окъ - Ma сащ 9,149. (58) 
1-1 р-О0 
и о). 
20 2 2 2 | 
в. = (С; + 81) (2A,/N,) (1 < PK) (59) 
апа 
Е асы И S (1 < 1 < РК) (60) 
1 лт — - 


иЦеге 5. апа с. are the data from the ith signal as in Chapter 


тув Ey C) denotes the expected value with respect to the PK 
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dimensional phase velocity vector with jointly Gaussian com- 
ponents as described in Chapter II. Utilizing Equation (15) 
from Chapter IV, some algebra and some trigonometric identi- 


ties, (58) can be written as 


K P-1 2 Р-1 Р-1 
СХ) ME т RR ‚+2 E 7, Е . R “ос 
V is] 2 28 рік Р 1 "20-21 рот ст 
тэ... (61) 
г раз = 


Ий (61), 


ф 


= principal value of (и М 


р.а акғі "pkti 


КТО Боркту: ЭЕ ea 
Arguments similar to those made in Chapter IV allow the 
high SNR maximum likelihood estimation equation for multiple 
stations to be written as 
KET- EPF I 


A(X) = E. expl X 2 7 
. i=] p=0 q=p+1 


2 


Е : 4 63 
Oi, p,q їр, ( ) 


In (63) e is the PK dimensional column vector with jointly 


Gaussian components just as defined in Chapter IV and 


Кркаі Как+а 
9. тарат 2 (64) 
эр» 2 (В ) 


рК-+1 t D 
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Пита БЕТА, ftormwmof "Equation (53) could be formed and the 
required expected value taken as in Chapter IV. The matrices 
involved will be quite large, however, since in general, four, 
five, or even more Omega station's signals are expected to be 
received. Since the solution to the estimation problem in- 
volves matrix inversions for each set of ET CAS. РИК Брине х — 
mom in an iterative process, the processing time required 
for the implementation of the optimum position estimator will 
be prohibitive. Although non-optimum, the most reasonable 
multiple station position estimation processor to implement 


is considered to be a pairwise station processor. 


Б/Н MULTIPLE ПОР POSITION ESTIMATION 

Multiple LOP estimation will require the estimation of 
the set of most likely LOP for the (P-1) different pairs of 
the P received stations. Consideration must be given to the 
spreading factor (defined in Chapter I) when choosing the 
Stations as pairs. The smallest spreading factors which can 
be achieved are desirable. 

Frenkel in [23] presents a maximum likelihood position 
estimation procedure which finds the optimum position esti- 


mate when given multiple LOP estimates. Although presented 


+ Abit hough some of the possible candidate components for 
X can be discarded before the matrix inversion is required 
by utilizing the “coarse check" procedure described in 
Appendix E, the processing time for the optimum multiple 
Station estimation will still be prohibitively large for a 
SAR application as shown in Appendix E. 
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as an ambiguity resolution technique, the procedure applies 
directly to the multiple station position estimation problem 
when the stations are treated pairwise as proposed. The co- 
variance matrix required for the position estimation proce- 
ноте сап вексај сем ated as in [23], 1f the variance of the 
phase estimates corresponding to the LOP estimates are known. 
In order to approximate what these variances will be, we 


recall Equations (25a, 25b) 





c, - 27m, (1 « i « K) (65) 


where d is the great circle distance of the receiver from 
the Omega transmitter generating the K signals. Since с. 

(1 < 1 < 2K) were assumed to be jointly Gaussian random 
variables with known covariance, the variance of each $i due 


to phase velocity uncertainty only is 


w.d\2 
с. 2 «( i ) а“. (66) 
ф. с Е; 


л. 





It has been shown as a result of [17, 21] that the noise dis- 
turbances to observed phase will in general be dominated by 
the phase disturbances due to phase velocity uncertainty. 

It is therefore a reasonable assumption to consider only 
phase velocity contributions in determining an idea of the 
variance of the 2K phases in the two station LOP estimation 


problem. 
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As a result of [17] and the experimental results re- 
ported on in Chapter VI, standard deviations of the phase 
velocity ratios (v/c) on the order of 3 x (ШЕ are expected 
for the four frequencies 10.2 kHz, 11 1/3 kHz, 13.6 kHz, 
and 10.88 kHz. Thus, Table VIII can be created, where d is 
the great circle distance from the station generating the 
signals to the receiver in nmi. If the standard deviations 
of the four frequencies given by Table VIII and the corre- 
lation coefficients given by Equation (49) in Chapter VI 
are used and the phases are considered to be equally weighted 


by the estimator, it can be shown that 


іде о 08 (67) 


~ 
~ 


JN да 


$ 


The value of O, in (67) corresponds to the standard devia- 


ф 
tion of the phase estimate for one station. The standard 
deviation corresponding to the LOP estimate resulting from 
a station pair would be the rms of the бү, associated with 
the two stations. The value of d required in (67) can be 
quite approximate and can be the great circle distance in 
nmi from the appropriate Omega station to the center of the 
presupposed area on the globe which was assumed given by an 
ambiguity resolution technique. The value of За вер Бу 
(67) is approximate, since the optimum estimator actually 
accounts for the estimated SNR of the signals instead of 


weighting them equally. However, (67) yields a result which 


is most probably adequate to enter the multiple LOP position 





TABLE VIII 


STANDARD DEVIATION OF PHASE DIFFERENCES 


i f, (kHz) ue d in nmi) 
125 10.2 E SS s 7 
2,6 11 1/3 11320 х 1073 
3,7 13.6 aora 
4,8 10.88 О 
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estimation procedure of [23]. Note that this б, only applies 
for four-frequency Omega using the frequencies specified in 
Table VIII. 

A final piece of information which must be provided to 
the multiple LOP position estimation procedure in [23] in 
order to obtain valid results is an effective phase velocity 
associated with the phase Q. This is required in order to 
do the mapping from phase to latitude/longitude on the earth. 
Although the phase velocities of the different signals have 
béen modeled as having known mean values with additive dis- 
turbances, it is expected that an average of the mean values 
of phase velocity will be sufficient to provide the necessary 
phase velocity for the mapping. These mean values will 
CES account for diurnal effects on the path being con- 


sidered. 
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This dissertation has reported the results of research 
efforts applied to the problem of achieving an optimum posi- 
tion estimator to be utilized for the GRAN worldwide SAR 
system. The principal new results which it contains are 
discussed below. Areas into which future research is needed 


are discussed as well, 


А. PRINCIPAL NEW RESULTS 

The dissertation has presented a new way to model VLF 
phase uncertainties. As presented in Chapter II, these uncer- 
tainties have been modeled as effective phase velocities 
which are jointly Gaussian random variables. The propagation 
effects on phase which can be reasonably well predicted such 
as diurnal effects are accounted for by specifying the mean 
values of the phase velocities. All the remaining effects 
which cannot be so readily predicted are modeled as contribut- 
ing to the uncertainty in phase through a randomness of effec- 
tive phase velocity. This model allows for very conveniently 
incorporating any correlations in the phase uncertainties, 
either interfrequency or interstation. This is accomplished 
by means of treating the effective phase velocities as being 
jointly Gaussian with known covariance. 

The maximum likelihood LOP estimator for an arbitrary 


number of Omega frequencies was derived for two stations in 
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Chapter IV. This estimator takes advantage of the model 

ERU S ued y discussed. In particular, the correlation between 
frequencies which is known to be significant and is presently 
actively being researched is simply incorporated into the 
estimator. If the interfrequency correlation which is ex- 
pected can be determined quite accurately, lane error occur- 
rence can be very significantly reduced. Uncertainties in 
phase or equivalently in phase velocity would still cause 
errors in LOP estimation but much fewer lane errors would 
occur if correlation between frequencies from an Omega sta- 
tion were high and predictable. 

As an integral part of the maximum likelihood estimator, 
doppler shift, signal amplitude, and noise density height 
estimation equations are derived in the appendices. The re- 
sults of four-frequency experimental data processed by the 
Signal amplitude and noise density estimators are presented 
in terms of SNR in one Hz in Appendix B. 

Monte Carlo simulations were carried out utilizing the 
maximum likelihood estimator derived in Chapter IV. These 
simulations were carried out for four-frequency Omega with 
fourth frequencies chosen primarily on the basis of ease of 
implementation at зе Omega stations, their separation from ` 
other knovwn or planned frequencies and the resulting unambi- 
guous ee size. As discussed in the previous chapters, it 
is nearly certain that for any of the known methods of 
ambiguity resolution which appear to have any hope of 


materializing for SAR systems, unambiguous lane size in 
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excess of the three-frequency, nominally 7/2 nmi lanes must 
be achieved. This fact prompted the simulations and it is 
now possible for the decision concerning the choice of a 
fourth frequency to be made in terms of lane error tradeoffs 
to achieve acceptable rane size. The fact that the simula- 
tions present results illustrating the extent to which lane 
errors can be expected as a function of SNR as well as phase 
velocity uncertainty and fourth frequency choice allow for 
observing the effect of SNR on lane errors. 

The simulation results prompt the conclusion that with 
adequate knowledge of the effective phase velocity variance 
and correlation and with reasonable signal-to-noise ratios, 
the unambiguous lane size can be increased beyond the present 
baseline nominal 72 nmi obtained with three frequency Omega. 
The controlling factors in the selection of a fourth frequency 
for Omega to be used for SAR purposes are (1) the capability 
of ambiguity resolution techniques to resolve the lane ambi- 
guity with restricted amounts of Omega data and (2) the 
accuracy with which the phase velocities of the Omege fre- 
quencies can be predicted. 

The maximum likelihood estimator of a line-of-position 
for the two-station Omega problem as derived in Chapter IV 
was applied to = и four-frequency Omega data. The 
results of this estimation (presented in Chapter VI) demon- 
strated that the maximum likelihood estimator was superior 
in terms of lane errors to the walk-up technique, [29], 


which is presently the only other method of processing the 
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experimental data. Of particular significance was the 20-25 
percent fewer lane errors experienced using the maximum like- 
lihood estimator for the second data period, night path con- 
ditions, over the walk-up technique employed in [29]. 

The experimental data processing using the maximum like- 
lihood estimator was done utilizing the best available ex- 
pected interfrequency effective phase velocity correlation 
coefficients (see Chapter VI) and various values of standard 
deviation for these phase velocities. The standard deviation 
values which yielded the best results in terms of minimizing 


2 


lane errors were on the order of 2 x 10 to 4 x pom as ex- 
pected from previous study by other authors [17]. The rela- 
tively large Ы of lane errors which were experienced 
for night and transition periods are suspected to be the 
result of higher-order-mode propagation interference on the 
Forestport Omega station path. The basis for this argument 
is presented in Chapter VI. A four-frequency experiment 
planned by the Naval Air Test Center, Patuxent River, Mary- 
land for late 1974 which is to relay the Omega data via 
satellite from various globe locations is expected to verify 
this argument.  Retransmitter (SARCOM) positions far enough 


from the Omega stations to insure small contribution from 


higher-order-modes of propagation are planned. 


В. FUTURE AREAS FOR INVESTIGATION 
Several areas into which future research efforts could 
be made have been observed as a result of this work, Since 


this is the first known attempt at modeling the VLF phase 
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velocities as random variables, prior efforts by other re- 
searchers have not analyzed the VLF data they have collected. 
in a fashion which readily yields the necessary parameters. 
The variances which can be "cd as a function of dis- 
tance for all the Omega frequencies need to be better deter- 
mined. 

There is presently some interest in investigating the 
correlation between phase uncertainties. Professor Pierce 
of Harvard has been examining the possibility of using the 
interfrequency correlation of Omega to advantage for some 
time. Lt is felt by this author that the determination off 
correlation, both interfrequency and interstation, ТЕ it can 
be determined to be significant even for some paths, is a 
most fruitful area for study. Interfrequency correlation 
could reduce lane error occurrence in Omega dramatically. 

In terms of the GRAN application, one of the areas which 
appears to have been neglected is the study of methods to 
reduce doppler shift of the signals received at the SARCEN. 
Proper design of the SARCOM and processing at the SARCEN 
could eliminate some of the sources of doppler shift and thus 
allow for an estimate of the velocity of the SARCOM to be 
made. In the present advanced OPLE concept, this was not 
considered possible Pelz For a SAR system, particularly for 
processing a distress incident at sea, SARCOM motion during 
and after transmission is a very real possibility. It will 
be much easier to locate a drifting life raft minutes or 
hours after notification if some velocity estimate is avail- 


able. 
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APPENDIX A 


DOPPLER SHIFT ESTIMATION FOR THE TWO=5TATION 


K-FREQUENCY CASE 


Ideally, if we knew the frequency of the received Omega 
signals at the SARCEN precisely (both wj and A. кота т. 
we could perform the necessary operations to obtain Equations 


(12a, 12b) in Chapter IV. They are 


10(j-1) + т, 


М 
с; = зА У: (Е) cos [ Co; Y A.)t]dt (A-la) 
210021) 
and 
" пао я 
5, = - б О ЖЕ К л оета (A-1b) 
: 10(j-1) 


The frequencies wj (1 < 1 SK) аге known quite well, however, 
because 2. A/R tone which is created in the SARCOM and 

sent with the Omega data to the SARCEN and because the frequen- 
cies at the Omega stations are well controlled. Therefore, 


we can form 


OO DET T. 


M 
ын ж== 
с. ne У; (Е) cos w,tdt (A-2a) 
: 10(3-1) 
апа 
" 10(j-1) * T, 
ї- 1 = 
5, ЕЭ у; (2) sin w.tdt (A-2b) 


10(j-1) 
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at the SARCEN from the received data, 
Tn order to form the likelihood function for A; 


К) Еве definitions 


A. т 


i 
Сі; = У; (Е) cos w,tdt . (A-3a) 
10:01:12 
and 
NOG укш 
Sij = У; (Е) sin Qitdt (A-3b) 


Шо > 


are made.  Recalling that 


У; (Е) = А. cos [ (w; + А ЗЕ - $1 + 07) + n; (t), (A-4) 


the probability density functions of Cis and 511 conditioned 
upon knowing 9. and Ф, can be found as follows. 
10(j-1) * T, 


! = -- 
С. А. сов [ (ш, 4 A.)t ф +0; ] соз w.tdt 
N9OCT-1) 


.10(j-1) + - 
t n; (t) COS w,tdt (A-5) 
10(3-1) 


can be written as 
- к T. 
(j-1) i 


А, 10 
г. Е = ' a 
с. = I. "c m + 0:)dt t E (A-6) 
1- 


where the second order frequency term is disregarded as not 


being in the passband of the system and 
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О «аз: л 


! = Е 
N13 | n, (t) cos wjtdt 5 (A-7) 
10(j-1) 


The integration in Equation (A-6) can be carried out yield- 


ing 


> 


нагі 
S ZA 


А {sin [10 A AT, - 9,3 P 
> MD A jer 5 0,1) + Nei (A-8) 


which can also be written as 








AI 
т FT : _ 

Ci j = z 9, cos ee + п. + Neij : (A-9) 

In (A-9) 

Sin (A.T.72) 
о, Е ЕАО. (А-10) 
быт 22 
ii 
na туш $i "m 9; з (А-11) 
y AT, 

and E 10 А. (3-1) + 5 : (A-12) 


Recalling from Chapter II that ni (t) is a sample function 
from a white Gaussian random process with two-sided spectral 
density N./2, and since Noid is formed as a linear function 


of n,(t), N'., will be a Gaussian random variable. N'.. will 
1 cij cij 


have mean value zero and variance 
10(j-1) + Т, 


а = Ми 5 ( 27. 4 2 
2 Т 2 OS DL Тоос 


10(1-1) 


NET 
1 
4 





E[N (A-13) 


cij 
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where the second order frequency term is again disregarded. 
Therefore, conditioned upon knowing фі and 01, с Ка 


Gaussian random variable with mean value 





1 = Í 5 u 
ЕС le и па (A-14) 
and variance 


МЕТ. 





VAR [ Ci; ] = (A-15) 


It can similarly be shown that 51. is a Gaussian random 


variable with mean value 





E[S!.] = - 


ij 5 sin LY i à + n.] (A-16) 


and variance 





VAR [ Sij ] = (A-17) 
Using the conditional probability density functions 
found above шин the fact that Ci; and 51. can readily be 
shown to be uncorrelated, and thus, independent since they 
are Gaussian random variables, the joint conditional density 
function of €. and 314 for the ith received signal can be 


written as 








M 
2 
1 1 = 
j=1 qu 
ANT 
2 ў 1 58 2 
ехр 1 NUS ИС, E П 2 


1 
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- Ai, 
1] 2 


2 


NT, 





2 
[S a, sin SEE F nod | ae (A= E 


This can be written as 








M 2 2А,0, 
, ' E De ee __- Хүү 
y er di 
AST. Q 
er Т 2 2 2 
= 3 Ни cue 1 + 1 
Sij Samy 73 + ту.) 2 | a 2 ЈИ 
МТ. А. ^а 2 ЗА, 0m M 
-ж ( Е коре шалны ] E n "EROGO Ecos Y - 
TN.T. 2N N. p ij ij 
01 2 j=1 
M 
1 。 ps : 1 1 : 2 
S sin Yas? sin n; Es COS Ya 4 + C sin 227 
М 
5, 2 (Сїй i зо an . (А-19) 
TI J : 


qe utilizing the definition of the likelihood function 
as given in Equation (9) of Chapter IV and again following 
the practice of not redefining the likelihood function when 
a monotonic function of it is taken because of the equiva- 


lence of the likelihood functions, 


ит А То“ дА а. M 
= шаг == A е 1 = 
A(A,) exp [ 2N. ] expí N, [cos ni mi cos Тір 


М 
1 . а ° 1 t A 
51. sin Үүл) sin n. ку cos 10 Ci sin Yij)]jf(ni)dni 


(A-20) 
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їп (А-20) £ (ni) refers to the probability density function 
of Пі: The random variables n, 0263 АҚЫ uniformly 
distributed on the interval (-T, T]. Therefore, Equation 


(15) in Chapter IV can be utilized to write (A-20) as 


ит А То” дА 0. Е. 
А(А.) exp [- e МЕ! 1, [==] (A-21) 
1 1 
where 
M M 2 
ВИС "с. сог у. "—_ ;  5:. віп Y..) 
1 ¿si ij = ij т 
M M 
Лаг” ЭТЭ ,Р uis o 1572 00:229 
іші 13 ij фар 13 íj 


The present OPLE ground station and the planned SARCEN 
call for processing the received signals at frequencies which 
result when the A/R tone is used for demodulating the folded 
or compacted spectrum., This will result in processing fre- 
quencies of less than 1 kHz. This method of processing re- 
sults in doppler shifts due to the local oscillator in the 
SARCOM which will not be removed. Thus, although because of 
motion of the SARCOM relative to the stations, the elements 
of the first K and last K sets of signals will be highly cor- 
related; this correlation is not taken advantage of and the 
doppler shifts are computed separately. Additionally, this 
method of processing, because of the reasonably unstable oscil- 
lator planned for the SARCOM units, will not allow velocity 


of the SARCOM to be predicted with any confidence [28]. 
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Reference [28] also contains a doppler estimation procedure 
which is similar to the results of this appendix although 
nd ddnde 

Evaluation of (A-21) over an expected range of doppler 
shifts in increments for each signal will yield estimates 
of the doppler shift of each of the 2K signals. These values 
of A. (1 < 1 <2 ЖК) сап Бе шеей tol fiind ithe valuepj of а. 
Ше 1 2 KOH and Үг (LA END ган ~ М) NSubsititutiorns 


of these values into the equations 


al: шигээ _ 
ij J ij cos Ya; a, ЗЕ sin 17 3 (A-23a) 
Ems Ci. Ein y + A. от (A-23b) 
ij Ын ij ij а, Е | 
M 
Ce == » C.. 3 (A-24a) 
a TET 33 
and 
M | 
5, = 2 S1. (A-24b) 
Базы 2 


ПОЩА i (I <i < 2K) yields thg desired doppler corrected 
data. The data elements С. and 5. can be shown, following 
the procedure in the first part of this appendix, to be 


Gaussian random variables with expected values 





= ti Ё ы 
Е (с. | E Фф.) 5 cos ( 9, + 9 1) (A-25a) 
and 
А МТ. 
ЕТ | 2017 2222: SUC DR cu (A=25b) 
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Their varlance will be 


N,MT, 
, $1) = 72. (A-26) 
aa, 


VAR (C, | 0,, ф,) = МАК (5, | 9 


i 1 
As the magnitude of the doppler shift of the ith signal in- 


Ф 2 Ф 
creases causing А.Т, 21 а, will approach zero as seen 


by (A-10). The variance of the data will get very large and 


the data will have an extremely small SNR. A standard зало 
curve will show that for AT, = t T, a,” = ор О. ТЕРІП 


other words, since Т. is nominally one sec. for all i, dop- 
pler shifts ofi magnitude < 0.5 Hzəcan be corrected with a 
deterioration in SNR of < A dB. A doppler shift of 0.5 Hz 
corresponds to a relative velocity of in excess of 10,000 
nmi/hr at the VLF frequencies, so with proper processing at 
the ground station utilizing the A/R tone, doppler shifts 
which are experienced in GRAN should be able to be estimated 


with little deterioration in SNR. 


140 


u 





APPENDIX B 


SIGNAL AMPLITUDE AND NOISE DENSITY ESTIMATION 


The maximum likelihood estimator derived in Chapter IV 
requires an estimate of the signal amplitude and noise spec- 
tral density height. It will be assumed that the data avail- 
able for this estimation is the received data corrected for 
doppler shift utilizing the estimation results of Appendix A. 
Therefore, the values of the Gaussian random variables ыг 
ап4 22 атг КЕШ ~ м) given by Equatson (A-23) will 
be the data. The value of K and M, as before, denote the 
nunber of frequencies and received Omega pulses respectively. 
The random variables EM and So. after correction for doppler 


shift will have conditional mean values and variance given by 











AT 
E(C,, | ВОИ) = — A = 6) (B-1a) 

271 
E(S,. | БА LED (B-1b) 

and 
NET, 
VAR(C,.|0,,6,) = VAR (s,,18,,9,)- о (B-1c) 
1 


In (B-1), А; is the amplitude, Т. the pulse duration, m the 
propagation phase delay, 9; the random phase shift and N./2 


the two-sided noise density height of the ith signal 
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(12-01: 2к Also ing B1) 


а = MMMM , ED 





where A. here is the estimated doppler shift of the ith sig- 
nal. The estimate of the value of A. and thus о. will be 
considered known for this derivation since it was assumed 
as the first processing step. 

The joint conditional probability density function of the 


data from the ith signal is 








2 
202 М 
ЕК ij EE Фі) таи x 
| а.“ М А 2 
exp {- ST (ға): (1; - 57 cos ni) + 





й 2 
Е 《Si 十 5 sin ni) 14) (8-3) 
Ще 
where 


This joint conditional density function can be written as 


2а : 
Y WM 
ij? Sij | DUREE. ‘т.т. 





1 1 
exp {- wo Es 4 2 - A, T, C, cos n, 
+ А, T, sin n, 17 (BS) 
where 
M 
NON 1 ler ро (В-ба) 
у= 4 j 
M 
C == 2 Can (B-6b) 
i sa 2 
and 
M 
зо беж о ооо eee (B-6c) 
1 ТЕ ij 


The definition of the likelihood function was given in 
Equation (9) of Chapter IV. Again, the practice of not re- 
defining the likelihood function when a monotonic function 
of it is taken will be followed because of the equivalence 
of the likelihood functions. The likelihood function for 


amplitude and noise density estimation can thus be written as 





T AA й 
БО он 
14 


1 
OT) т foo = 
N, 
1 
2а 2А 


i i 
exp | q 


| (с. cos n, - S 
i 


; Sin "121. 


£ (ni) dn, Я (B-7) 
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£(n,) in (B-7) is the probability density function of пі» 

In Chapter 11 9, ПЕ Е УУ шеге shown to be uniformly dis- 
tributed on (-T, T], so all the n; Vill also have this den- 
sity function. This and Equation (15) in Chapter IV allow 


the likelihood equation to be written as 





1 20.“ A, Mr,” 
N ie 
i 
20. ^A. B. 
(=== (B-8) 
O N 
i 
where 

2 DELL 2 
В; (с. + 5. ) : (8-9) 


rio) in (B-8) is the zero order modified Bessel function of 
thè irst kind, 
For estimation of А, the likelihood function is written 


as 


дал 1 202172 
1 1 1 1 


) I (———. (B-10) 


А(А,) = exp (- 2N, o N 


Taking the partial derivative of (B-10) with respect to А. 


and setting the result equal to zero yields 


2а А.В. 
T ) 
= 28. 1 N, 
А. 二 一 一 А (B-11) 
= нг а^ Б. 
T ( 1 ==; 
O М. 
і 
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where А, is the estimate of A, and 1,0) is dhes irst ‘order 


i 
modified Bessel function of the first kind. 


For estimation of N the likelihood function is written 





12 
а5 
2а 2 А Тит? 
AQ) = expl QA (Е, + o . 
N 
i 
22, ^4; | | 
1 е Е (B-12) 
O N, 


Taking the partial derivative of (B-12) with respect to м. 


and setting the result equal to zero yields 


2a ГА В. 
2 2а ТЕ. EA cT 2а Ep c E N, 
Te ee O a 1 
i шл р: 
І ( 1 і, 
O N. 
1 


where N is the estimate of Ni. 

Equations (B-11) and (B-13) do not have the values to be 
estimated isolated and further, they must be jointly esti- 
mated since each estimate is 111 upon the result of 
the other. A convenient method of jointly solving the equa- 
tions isi to isolate? thel а Ето Ii ATE (°) in each equation 


and equate them. This yields 


^ zu = A, i 19 
М, = ур T 2 : (B-14) 
i 
WAS 





Substitution of this value of N, into (B-11) yields 


4A.B.MT. 
ТО 1 
їі ‹ 4 2. 2) 
^ 28, 46.-МА, T, 
К ут ЖЕ ажа å (8-15) 
E i 4A.B.MT. 
1 ( Hi E 1 ) 
O 2 2 


A EA өш 


Equation (B-15) can be iteratively solved for А. now and us- 
ing the resultant value of A, in (B-14) yields an estimated 
^ 
value for N.» 
In order to discover a method convenient to solve (B-15) 
iteratively, the following was done. For relatively small 


NIE in other words high SNR, from (B-14) 


28146. А то” 
вит ==. (B-16) 
1 
От МА тт.” 
в 241, (B-17) 


which implies the argument of the Bessel functions in (B-15) 


are large. For large arguments, ТС) and I) are approxi- 
mately equal [27]; therefore, for high SNR 


28, 


Тм ° 
i 


> ~ 
шог 
1 





(8-18) 


‚> 
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Epis value, is#rormed®@fron the data Гог еасп 1 (1 < i < 2K) 


and definedhas 2n initial condition for an itëērative procēss 





A 28, 
т эг 
1 
the (N-tl)st iteration vialue*for А. can be shown to be 
DNE м?т, 
ы E 2) 
А rw Үг" сн 
i, N+1 i,o 23 2 м?т 2. 
i,N 1. 
M ( a 2 2) 
BT ü д N m 


A computer program was written to do the iterative pro- 
cessing in (B-20). By substituting the resulting estimate 
of А; into (B-14) an estimate of М. мас obtained. The flow 
chart of the computer program is contained in Appendix E. 
From these estimates, the signal-to-noise ratios in one Hz 
which the Omega signals have at the front end of the SARCOM 


receiver can be estimated from 





ANN T: (B-21) 


The results of evaluating (B-21) with the four-frequency 
data obtained from the experiment discussed in Chapter VI 
are presented in Figures Bl - B24. Because of the signifi- 


cant diurnal effects on phase velocity, as discussed in 
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Chapter II, the data has been divided into three major cate- 
gories. Both transmission paths in daylight is referred to 

as "Day" condition, both paths dark is referred to as "Night" 
condition and partially lighted paths is referred to as 
"Transition" соса ценена The SNR estimation results of the 
first five-day data period were so little different from 

the second five-day period that the two periods have been 
combined. The power radiated by the Forestport Omega station 
for the data collection periods was varied over quite a wide 
ránge as weather conditions varied.  Radiated powers on the 
order of 100-200 watts for the lower three frequencies and 
200-400 S for 13.6 kHz were typical. The radiated power 
from the Trinidad station was on the order of 1 kW for all 


frequencies. The finally implemented Omega stations are 


expected to radiate 10 kW at all frequencies consistently. 
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APPENDIX C 


LOP ESTIMATION FOR KNOWN PHASE VELOCITIES 


In developing the model for the estimation problem to 
be solved, the first approach which was used was to assume the 
phase velocities were well known. The values of phase veloc- 
ities which were intended for use were those given in Table 
ТТ, the first mode propagation values obtained from Professor 
I TA. Pierce of Harvard University. The derivation of the 
likelihood function for this model is carried out in this 
appendix. Simulations with several possible fourth frequen- 
cies were carried out and they are presented here to serve as 
an upper bound for the performance which can be expected from 
the maximum likelihood estimator when the uncertainties in 
the phase velocities are accounted for in the derivation. 

The difficulty with assuming well known phase velocities 
is that as the signal-to-noise ratio gets reasonably large, 
the estimator places great Mont idence in the data received. 
Since the additive noise is the only factor disturbing the 
data in the known phase velocity model, as SNR increases, 
very good data would be expected. Literature search and some 
preliminary experimental data processing revealed that the 
assumption of known phase velocities is not justifiable in 
a model and thus the model described in Chapter II was created. 
Because of its value in bounding the performance of the esti- 


mator and in obtaining a feel for the problem, the known 
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phase velocity estimator is derived herein. 
As was pointed out in Chapter IV, Equation (16) can be 
used as a starting point for the known phase velocity esti- 
mator derivation since no approximations or operations in- 
volving the phase velocities' randomness have been performed. 
Therefore we have as a likelihood function 
2 2 


К 
x) t І. LER; О. 28, К „1 Соз Y] 


1/2, (с-1) 


where the variables in (C-1) are defined as in Chapter IV. 
For high SNR, the same arguments as applied in Chapter IV 
can be applied. The likelihood equation can thus be written 
as in Chapter IV, Equation (30) as 
S 2 
A(x) = exp [ 2. Q, 9,7] (с-2) 
121 
since no application of the randomness of the phase velocities 
was made to that point. In (C-2) Qi and Ui are defined as in 
Chapter IV, Equations (17) апа (31), The log-likelihood func- 
tion, an equivalent likelihood function to maximize, is 
convenient to form here, particularly for computer evaluation 


to avoid overflow problems. Therefore, 
NX) = | Шы 2 (C-3) 


Monte Carlo simulation for several fourth frequencies which 


could easily be added to the Omega format was carried out 
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utilizing the estimation Equation (C-3). The results of the 
simulation in terms of the probability of lane error occur- 
rence are presented in Figure Cl. The definition of lane 
error here is as defined in Chapters IV & V. 

Two general once can be drawn from the simulation 
assuming known phase velocities. First, since there is a 
degradation in probability of lane error on the order of .1 
to .2 in going from unambiguous lane size of a nominal 500 
nmi to approximately 1000 nmi, the latter size can serve as 
a soft upper bound in selecting fourth frequencies to do 
simulations when jointly Gaussian distributed phase veloci- 
ties are assumed. Second, a fourth frequency of 12.364 kHz 
which results in a nominal 785 nmi unambiguous lane appears 
to be worse in terms of lane error probability than either 
10.462 kHz or 10.737 kHz which result in nominal 927 nmi and 
1356 nmi unambiguous lane sizes respectively. This implies 
that frequencies chosen from the lower part of the 10-14 kHz 
spectrum can be expected to yield better estimation results 
than those chosen from the upper part of the spectrum in 
terms of lane error probability. This can be explained by 
the fact that the ith frequency being higher yields adjacent 
maxima of the ith component of the likelihood function which 
are closer together and thus the chances of small perturba- 
tions in the phase of the ith signal caused by noise will 
more readily cause the peak of the likelihood function asso- 
ciated with the incorrect LOP to be larger than the correct 


peak. 
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APPENDIX D 


EXPECTED VALUE OF A SPECIAL FUNCTION 


Chapter IV derived the likelihood equation given by Equa- 


tion (32) as 


AGO = Е ехр [е Ае + Bue + D] . (D-1) 


In (D-1) Er) denotes the expected value with respect to € 
where € is a column vector with 2K elements jointly Gaussian 
Ek buted with zero mean values and covariance matrix Ko: 

A is a 2K square matrix, B is a column vector with 2K elements, 
and D is a scalar. The expected value in (D-1) can be written 


as 


A(x) Eu xr та / 45520 / expí- [ЕТ Ko os 
Т К 


оос 2e Ae Be 461,...... de (D-2) 


2K 





where E denotes the determinant. The argument of the expo- 


nential in the integral defined as G can be written as 


G = ЕТК ас Е - 2є ЇАє - 2BTe = eT (k 71 - БАЈЕ-В'Е-Е ТВ, (D-3) 


which can then be written as 





G = о MANE ик“ -2A)E -E OK 7} AD) Y (D=4) 
where 
eee (5 7-24)” (D-5a) 
and thus 
U = FAS M NN (D-5b) 
Equation (D-4) can now be written as 
С = (є - E ex 7 - 2A) (e - 0)- шилэн - 2A) U (D-6) 


by completing the square. 


These manipulations allow (D-2) to be written as 


ў ехр [р + iv (k 7 ЕКО A eau) 


(х) ч ------------------2------. 
(21735 | x, | 1/2 


[ 522: | expí- 5 [(є-и) T ско -2A) (e-U) ] 3 de, “ЧЕК : 


Е Е 


1 2K 


US, 


The integral is over all values of the jointly Gaussian random 


vector € and thus equals 


1 


(лук | (к - р а 
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The desired expected value is then 


| (NES 2A) 1 |842 en Е 
A(x) = —— exp [D+ 20 (Ko -2A)U] .(D-8) 
K 
E 


Matrix manipulation of this equation allows the final simpler 


result 


exp [D + ві Киа MARD B] 
A(x) = a a a ЕЕЕ (D-9) 
BE кк 


to be written. 
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APPENDIX E 


COMPUTER IMPLEMENTATION OF TWO-STATION LOP ESTIMATION 


In order to solve the estimation equation given by Equa- 
tion (38) in Chapter IV for the two-station, four-frequency 
experimental data discussed in Chapter VI, a computer pro- 
gram was written. A flow chart of that program is contained 
in this appendix. On the IBM 360-67, the program occupied 
34K BYTES of core. The time required to process one "run" 
of data was on the order of 15 seconds. A run, as defined 
previously, consisted of nominally three minutes of Omega 
data. The range of values over which x was allowed to vary, 
as discussed in Chapter VI was twice that which would be re- 
quired in an operational system, 

A "coarse" processing which was performed to save com- 
puter time was the following. The values of Y, Qd < i < 4) 
as given by Equation (17) in Chapter IV were examined with 
the value of x being tried in the iterative procedure before 
entering the actual el ue manipulation portion of the pro- 
gram. The magnitude of each value of m was required to be 
less than 1.5 rad. before the entire estimation process was 
carried out. This is equivalent to saying that noise and 
phase velocity uncertainties will not disturb the phases by 
more than + 1.5 rad. for the Dallas, Texas receiver location. 

The value which is used for this "coarse" processing can 


be determined for the ith frequency by examining the amount 
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of phase shift which two standard deviations of phase veloc- 
ity disturbance from each station in opposite senses will 
cause. To determine the effect of a phase velocity distur- 


bance on phase the equation 


- мэ х 
tk A A ey 
К-т E 


СЕЗ) 


is written. The values of d are the propagation phase delays 
for the ith frequency signals, o. is the ith angular frequency, 
VQ. are the phase velocities, m. is the number of whole cycles 
in the phase difference and (L-x) and x are defined as in 


Chapters IV and VI. To see the maximum effect of perturbing 


phase velocities, let Av be added to Vers and subtracted from 
У and for simplicity, assume У аи e інн У + Thus, 
2, Pros x 
Ze, ac а Zu 
ER L-x х 
= — (------ ~- —) - 2Tn. (E-2) 
У 1457 у AY 1 
у у 


where Ab, is the perturbation in the phase difference caused 


by Av. This can be written as 





- - А = Трах) E аен у 2 
р A х mt APRI aan. 
ші wb Ау 
m (L-2x) - атш оса) 95550), 


since (E, will be small. Or, 
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At Av | E i | Av | 
~ . 








149, | Е (E-4) 
For L = 3547.15, which is the Forestport-Dallas-Trinidad 
Av 
distance, о, = ZO OZ Гапа | 一 | = 8 x 10 3; the magni- 
С 
tude of the maximum expected phase perturbation is found to 


be 


ои on 


~ 2 = = 
[Ago] = 1618757 1.5 жый . (Е-5) 


Since the | 49 | for lower frequencies will be smaller, this 
was considered a reasonable upper bound to use for all fre- 
quencies of the four-frequency experimental data. 

The amount of savings in terms of the number of values 
of x which must be considered when a coarse search such as 
this is used can be determined as follows. Consider the ith 
frequency phase difference as x is varied through a lane. 
Only one half the values of x vill meet the coarse criterion 
if ly, | is required to be less than 5. Since there are 
four frequencies, only 1/16 of the possible x will meet the 
Coarse critierionsss As discussed im Chapter УГ, + 350 nmi 
(nominally) were used in processing the four-frequency experi- 
mental data with increments of x = 0.5 nmi. Thus, only about 


90 values of x vice about 1400 had to be processed through 


the entire estimation procedure. 


170 





For a generally implemented system such as GRAN, it would 
not be unusual to have values for L = 7000 nmi or greater, so 
it is not expected that this processing advantage can be 
generally implemented. Because of the significant time 
savings it affords, it may well be worth using a variable 
value which is determined after a rough estimate of L is 
known. In the final GRAN system, if 10.88 kHz is added as 
a fourth frequency, 360 nmi unambiguous lanes will result. 

If Ax of 0.5 nmi are used, 720 values of x must be processed. 
ТЕ по eoarse jesitimate: procedure is used, usingithe ^ 15- 
second processing time for twice the number of values of x 

as a ra} about two minutes will be required to estimate 

an LOP for a station pair using a program similar to the 

one in this appendix. 

Preferable methods of search for the values of x which 
maximize the likelihood function may exist, such as [30]. 
Methods other than the exhaustive search method discussed 
in previous chapters of this dissertation and implemented 
in the following computer program were not considered. The 
exhaustive search technique was chosen because of its sim- 


plicity and ease of implementation. 
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v(I)=FR2-CL(I 
+(1-FR2)CD(I) 
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2R(I)*SQ(I) 


| Write 
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Q(I) = 
-R(I)-R(I+4) 


| 2R(1)+2R(144) 










174 





— Ро 





I=1 
! 
[ee | вм FAS 
| 
1-1 





175 








на 
| 
Е 





: 
| 
нн 
2 
À 
Ен 
= < И = == 
~ | ЈУ Ы 
М 
ү 


D-D+W° (T)Q(T) le а 


2-20 : WI) 
(ш«(1)/с) 


В(1)--Х:2 


BI.) = 
І-Х):2 


А(1,1)-Х%22 


A(I+4 1+4) = 
(L-X)2.22 


E 
а, 


176 





A(I,I+4) = 
-Х(1-Х)-22 


А(144,1) = 
A(I,I+4) 


{= 


LU 


* The operations in these boxes 
are matrix operations. The sym- 
bol |° [| denotes the determinant, 
i dengtes the identity matrix 

and B' denotes the transpose of Y 


th tor B. | 
e vec r 1-(|1-2А5,| 7 
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The parameters to be read in are the scalars XS, XF, AX, 
L, FRl, FR2, and IDIM and the matrix к. 


Tef invariants to bef set up arelthe scalars T, RNG and c 
(vel or тег in vacuun end the? vectors .T, ш, со, ср 
and RNG. 


The (I,J)th data element consists of  ec(I,J) and es(I,J) 
where 


10(1-1) +Т (Т) 


ес(Т,1) = / у(1) COS[w(I)t]dt апа 
10(J-1) 
10(J-1)+T(I) 

ec(I,J) = f y(I) SIN[w(I)t]dt 
10(J-1) 


where y(I) is the Ith component of the received signal 
vector. 


Subroutine NEST must be provided the vectors R, Б апа Т 
and returns in R the estimates of the signal amplitude 
to noise density ratios. 


The vector SWC must be provided as propagation correction 
factors for each frequency in the form of station two 
minus station one corrections., When propagation correc- 
tion factors are used, all entries in the mean value 
vectors V must be set equal to the charting phase velocity 
which is used in their prediction (1.00261*c). When no 
propagation correction factors are used, the vector SWC 
must have all zero entries. 


Subroutine ORDER must be given the values of the scalars 
U, X and IDIM and it returns the requested number (IDIM) 
of most likely estimates of X in the vector XHOLD. The 
corresponding likelihood function values in the vector 
HOLD are returned also. The components of these vectors 
are in descending order of likelihood. 


The vector HOLD at this point contains the probability 
ratios. 
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DEFINITIONS 


Х5- starting value for X in the iteration 
XF- final value for x in the iteration 
AX- size of incremental steps in X to be used 
L- fixed value for L | 
Kon 8-by-8 covariance matrix of € 
CL- 4 element vector for mean values of phase velocity 


in lighted path 


CD- 4 element vector for mean values of phase velocity 
in dark path 


T- 8 element vector for Omega pulse durations 
FR1- fraction ог pathi Imin light 

FR2- fraction of path 2 in Wight 

0- 4 element vector of radian frequencies 


Y- 8 element vector of mean values of phase velocity 
to be computed by the program 


SWC- 4 element vector of propagation correction factors 
in form of station two minus station one corrections 
for each of four frequencies (to be set equal to 
zero if no corrections are to be used) 


RNG- a predetermined constant in radians equal to the 
maximum deviation of the phase of any signal from 
the phase resulting from mean values of phase 
velocity. These deviations are due to phase velocity 
uncertainties. (1.5 rad was used in the experiment 
data processing.) 


IDIM- The number of most likely estimates of X to be 
retained 


HOLD- vector of length IDIM of values of the likelihood 
function to be retained 


XHOLD- vector of length IDIM of values of X corresponding 
to the values of HOLD 


u- 8 element vector of phase angles computed by the 
program 


181 Е" 
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Subroutine NEST computes the joint maximum likelihood 
estimates of the signal amplitudes and the one-sided 
noise spectral density heights from the equations 


10 _ SOG a 


9T(I) 2 


“(1 “ши . 1|12А(1) 550 (1) /Х(1))1 


9-T(1) 1 [2A(1)*SQ(I)/N(T)] | 


N(1) 


In these equations the vectors are as defined in the 
subroutine NEST or the MAIN estimation program. I, and 

I are modified Bessel functions of the first kind of 
order one and zero respectively. The ratio of these 
Bessel functions is computed by using their polynomial 
approximations found in Handbook of Mathematical Functions 
by M. Abramowitz and I. A. Stegun, Dover Publications. 


The constants required for this subroutine are: 


Cl = 2.50662828 Dl = 25339353553 
C2 = .39894228 D2 = 1.00000000 
С3 ® ‚053416616 D3 = 175781160 
C4 = „01118812 D4 = -1.02993610 
C5 = -.00161279 D5 = . 90497756 
C6  « 01920037 D6 = -.84749981 
С? = -.04017321 D/ = .80494379 
C8 = .04837180 D8 = -,76696533 
C9 = -.02678420 

ClO = .00532070 


It has been found experimentally that 10 iterations are 
adequate for reasonable convergence of the iterative 
process for data from the OPLE ground station in SNR 
environments which can be expected from Omega. 


The vectors A, ZN, SNR and R written out here are the 
estimates of signal amplitude, two-sided noise density 
height, signal to noise density ratio and signal ampli- 
tude to noise density ratio respectively. 
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* HOLD, XHOLD;X,U,IDIM 








187 





HOLD(I-1) 
= HOLD(I) 


XHOLD(I-1) 
- XHOLD(I) 





| Subroutine ORDER receives the scalars X, U (likelihood function value 


for X) and IDIM. ORDER places U and X into the vectors HOLD and XHOLD 
of dimension IDIM in their appropriate place. That is, the largest 
value of U received thus far in the "run" occupies the first component 
of HOLD, the second largest U occupies the second component, etc. 

The components of XHOLD are the values of X which correspond to these 
values of U, 


2 This operation insures that only one value of U and X will be retained 
for each peak of the likelihood function. 
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